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Preface 



The main purpose of this work is to study the perturbed integrable systems in the semiclassical 
regime |33|. The examples of such systems include a circular or rectangular billiard with a 
distorted boundary a rectangular billiard with a weak magnetic field, coupled quartic oscil- 
lators — any system that is a "small" perturbation of an integrable system. In a classical 
integrable system even an infinitesimal perturbation produces a qualitatively big effect — it 
changes the topology of the phase space near the periodic orbits. Whether this change will 
affect the quantum states depends on the energy of the particle. 

As will be explained in due course, these systems possess the special spatially localized 
states related to the classical periodic orbits. Localization may take place even for relatively 
small perturbations. For example, in the case of a distorted boundary even a perturbation 
smaller than the wavelength may lead to a strong effect. 

We use the Bogomolny T-operator method throughout this work. It provides the semi- 
classical description for the Poincare surface of section. With this approach we reduce the 
dimensionality of the problem by one. In particular, for the two-dimensional systems one can 
derive a one-dimensional Schrodinger equation on the surface of section. It is often easy to 
find the qualitative behavior of its solutions and thus predict the localized states even without 
doing extensive calculations. In some cases the T-operator method produces the results similar 
to the Born-Oppenheimer approximation that we also consider. 

Here is a brief description of the work. Chapter 1 is of an introductory character. We 
review the T-operator method in general for the reader's convenience. We also demonstrate its 
connection with the boundary integral method for billiards. The Maslov phases are discussed 
as well. 

In Chapter 2 the perturbation theory is derived and discussed from various aspects. Possible 
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experiments are suggested. Although the theory is formulated for the circular billiard with a 
perturbed boundary, the general character of the derivation and results is emphasized. In the 
following chapters we apply the method to other systems without detailed explanations. 

In the first part of Chapter 3 we adapt the theory to the rectangular billiard with a 
perturbed boundary. We consider a tilted square as an example. In the rest of the chapter 
we study several non-perturbative systems where the perturbation theory still can be used for 
special classes of states, perhaps after some modification. Among such states are the bouncing 
ball and the whispering gallery modes. We also show how the whispering gallery mode near 
the boundary with a point of zero curvature can be described in terms of scattering. 

In Chapter 4 we consider an example where the perturbation is not a distorted boundary, 
but a magnetic field. We analyze the square billiard in uniform magnetic field and with an 
off-center Aharonov-Bohm flux line. Experimental possibilities are also discussed. 

The Born-Oppenheimer approximation is a subject of Chapter 5. We revisit some of the 
examples of earlier chapters and analyze them with this method. We compare this approach 
with the T-operator perturbation theory. 

In Chapter 6 we derive the trace formulas for chaotic, integrable, and almost integrable 
systems starting with the T-operator. We consider an example of coupled quartic oscillators. 
Finally, in Chapter 7 we summarize the results of this work. 
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Chapter 1 



Introductory chapter: Elements of the 



semiclassical theory 



We begin with a review of the semiclassical methods that provide the framework for the 
perturbation theory developed in the subsequent chapters. The starting point for most of our 
calculations will be the Bogomolny T-operator equation, which is essentially a semiclassical 
Green's function method adapted to the Poincare surface of section. In billiards the procedure 
is related to the boundary integral method. 

1.1 Stationary phase approximation 

Almost any semiclassical theory takes advantage of the stationary phase approximation (S$>), 
which provides the most direct connection with the underlying classical behavior. Since it 
is used consistently throughout this work it seems necessary to remind the reader of the 
application of this method. The problem that often arises is to evaluate an integral of the 
form 



asymptotically for h — > 0. It is argued that under certain conditions the main contribution 
comes from the neighborhoods of the points where the derivative of S(x) vanishes, the so called 
stationary points. Near such point x st we can expand 




(1.1) 



S(x) = S(x st ) + \s"{x st ){x - x st f + O [(. 



X - Zst) 3 ] 



(1.2) 



1 



and the integral is approximately equal to 



The square root is analytically continued from the real positive numbers. Equation (h3) is a 
good approximation to the integral if the third order terms that we left out are small within 



the region \x — x st \ < ^Jhj \S"{x & t)\ that contributes the most to the integral. This translates 
into the condition 

|5"(x st )| 3/2 »ri 1 /2|5»'( a;st )| (i.4) 
which is normally satisfied for small h unless |S"'(a;st)| is unusually small. 



1.2 Boundary integral method 

The perturbation theory that is the subject of this work is based on the semiclassical surface 
of section transfer operator (or T -operator) method developed in generality by Bogomolny 
Jil| and discussed below. This method is especially convenient in billiards where the surface 
of section can be associated with the boundary. It is instructive to see how the T-operator 
formulation follows explicitly from a version of the boundary integral method JT^, |33| designed 
to solve the Helmholtz equation in billiards. 

The wavefunction ^f(r) for a particle of mass m and energy E moving freely inside a two- 
dimensional domain B with the impenetrable boundary dB is determined by the Schrodinger 
equation 

(V 2 + k 2 ) = (1.5) 

with the Dirichlet boundary conditions 

*(r') \ dB = 0. (1.6) 

Here the wavenumber k = y / 2mE/H and V' 2 is the two-dimensional Laplacian acting on r' . 
The free-space Green function Go(r' , r) (that does not satisfy the boundary conditions on dB) 
is defined by the equation 

(V' 2 + fc 2 )G (/,r) =S(r'-r) (1.7) 
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which is solved by the Hankel function J41J 



G (r',r) = - t l H^ (k\r'-r\). 



(1.8) 



Multiplying Eq. (1.5) by — Go(r',r) and Eq. (1.7) by 5 , (r'), adding them together and inte- 
grating over the billiard's area we arrive to the integral equation 



*(r) = / dV [*(r')V ,2 G (/,r) - G (r', r)V' 2 *(r')] . 
Jb 

With the help of Green's theorem we transform it to the boundary integral 



*(r) 



dq' 



(1.9) 



(1.10) 



dn' dn' 

where q' is the coordinate along the boundary and n' is the normal at point q' directed inside 
the boundary. The second term vanishes due to the boundary conditions. Taking point r = q 
on the boundary and introducing the new function 

d$(q) 



Kq) = 



dn 



(1.11) 



we can write the integral equation which is equivalent to the original Helmholtz equation ( |l.5| ), 

,dG (q',q) 



fi(q) = / dq' 



di 



W) 



Its kernel has a ^-function type singularity at q' — q. Indeed, 

dG (q',q) ik dL(q,q') 



with 



dL(q,q') 

dn 



\P -n\, q' ^q 

i. q' = Q 



(1.12) 



(1.13) 



(1.14) 



where L(q,q') is the length of the classical orbit (chord) going from point q' to q and p is a 



unit vector along this orbit (or a unit momentum) (Fig. 1.1). Using the asymptotic expansion 



for the derivative of Hankel's function for small argument 1 41 we find 



dG (q',q) 
dn 



27rL(g, 9 ') 

oo, q' =q 



(1.15) 



The — 1/4tt is written under the assumption of a convex boundary with non- vanishing curva- 
ture, but this assumption is not essential for extracting the (5-singularity. So, we may assume 
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Figure 1.1: The classical orbit from point q' to point q on the boundary. The dashed line 
shows the orbit incident at q' . 



that the boundary is locally a straight line. Then, for any q and q', dL(q,q') / dn =lim 



n I ' \j (q — q 1 ) + n 2 and the singular part of the kernel 
'dGo(q\qY 



dn 



sing 



= lira 

TWO 



2tt 



(q-q'Y+n 2 



The integral equation ( 1.12| ) may now be rewritten as 
with the regular kernel 



dq'K(q,q') f x(q' 
fdG (q',q) 



(1.16) 



(1.17) 



(1.18) 



\ dn 

Thus the problem is reduced to the one-dimensional integral equation ( 1.17 ) for u(q). This 
function contains most of the information about the state and det(l— K) = is the quantization 
condition. [The second quantization condition comes from the periodicity requirements on 



fJ-(q)-] If needed, the wavefunction \l/(r) can be subsequently determined from Eq. (1.10). 

Our goal is to develop a semiclassical approximation to this scheme. When kL 3> 1 we may 
use the large argument asymptotics for the Hankel function pdj| and the kernel now becomes 



K{q,q') 



2TrL(q,q 



7j \p-n\e 



ikL(q,q ) — z7r/4 



1_ d 2 S{q,q l ) \p-n\ iS[q , q ')-^/ A 



V 2-7T dqdq' \p' ■ n'\ 
We introduced the reduced action (measured in units of fi) 



(1.19) 



S(q,q') = kL( q ,q') 



(1.20) 
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for the classical orbit q' i— > q. Clearly, K(q, q') is not symmetric. The semiclassical approxi- 
mation allows to symmetrize it. First, notice that, in spite of the notation, both \p ■ n\ and 
\p' ■ n'\ are the functions of two coordinates, q and q' . The essential part of the semiclassical 
approach, however, is to evaluate all the integrals containing the fast changing action at the 
exponent in the stationary phase approximation. It will be clear from the following chapters 
(and is well known) that this approximation amounts to extracting those q and q' that are 
connected by a classical orbit. Thus we may consider, say, \p ■ n\ to be a function of q only 
(as well as the two fixed parameters related to the initial conditions, one of which is the total 
energy). Similarly, \p' ■ n'\ will be a function of q'. 
Define a function 

m - -0= r (1.21) 



Then the semiclassical analog to Eq. (1.17) is 



#?)= / dq'T(q,q')i;(q') (1.22) 

JOB 



where the T-operator 



Unlike K(q,q'), the T-operator is symmetric and unitary (semiclassically). See Sec. 1.5 for 



more details. Clearly, Eq. ( 1.22 ) has a solution if det(l — T) = 0. The original wavefunction 



^f(r) can be found from ip{q) (see Sec. 1.4). We emphasize once again that the stationary 
phase approximation is expected when dealing with the T-operator. In the next section we 
discuss the relationship between the T-operator and the classical surface of section map. 



1.3 Surface of section map 

In the previous section we reduced the two-dimensional problem in a billiard to a one-dimen- 
sional problem on its boundary. The state of the system can be equivalently described by 
either 2D-wavefunction \&(r) or ID- wavefunction ip(q). In the classical language, we reduced 
the continuous motion in the four-dimensional phase space of the billiard to a map in the 
two-dimensional phase space of the boundary It is an example of a surface of section map. 
Again, this map contains all the information about the continuous motion. 
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Formally, the map (p', q r ) i — > (p, q) is implicitly given by the equations 

dS(q,qi) , dS{q,q>) 
P = ~ > P = 



dq' 



(1.24) 



The generalized momentum (in the units of ft) p is the projection of the total momentum onto 
the boundary. To connect the classical and quantum descriptions consider a version of Eq. 



(Il-22|) 



i>(q) = J dq'T(q,q')i)'(q') (1.25) 

describing the propagation of a wavepacket. We write the wavefunction ip'{q') on a locally 
defined Lagrangian manifold p'{q') in the semiclassical form (Sec. 6.3.C in Ref. 173]) 



ld*S(P,q>) JJp , ww = IfdQ 



Jf p'(q')dq' 



(1.26) 



V dPdq' y \dq'y F 

where S(P,q') is the Legendre transform of S(q,q') to the canonical coordinates (P,Q) such 



that P is a local integral of motion. We substitute it in Eq. (1.25) and integrate by the 5$. 
The stationary point q' at is determined by the second of the Eqs. (1.24). In the neighborhood 
of point q we define the function p(q) by 



q pqL(q) 
p(q")dq" = S[q,q' st (q)}+ / v\4'W ■ 



(1.27) 



This definition satisfies the first of the Eqs. (1.24). Thus the integral 



\ 



d 2 S 




dqdq' 


{oq')p 


d 2 S 


I d P' 




r dq' 



Ppif P{q)dq 



dQ 



J J P(q)dq 



ld 2 S{q, P) 



e ijp(q)dq-in = (1. 2 8) 



p ywy / p W OqdP 

(To prove the second equality differentiate the equation dS(q, q')/dq' +p'(q') = 0. The Maslov 
phase — 7r results from the reflection at the billiard's boundary.) This shows that if a manifold 
p' {q') evolves into the manifold p(q) the T -operator provides the evolution of the wavefunction 
defined on this manifold. 



1.4 T-operator: the general case 

The phase space of billiard boundary is an example of a Poincare surface of section (PSS). 
The above discussion makes it clear now why the T-operator is called the surface of section 
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transfer operator. In general, PSS is a (2N — 2)-dimensional manifold in the phase space 



of an TV-dimensional system crossed by all classical trajectories. Bogomolny |19 derived the 
expression for the T-operator in the case when PSS is an (N — l)-dimensional manifold in 
the coordinate space of the system together with its conjugate momenta. He assumed the 
Hamiltonian has the form 

(1.29) 



H(p,r) = -f + V(r) 
but presumably his result is more general. The T-operator is now 



T(q,q';E) = 



1 



,{N-l)/2 



E 



rd 2 S(q,q';E) 



det 



1/2 



JS(q,q';E)-ifu 



(1.30) 



The sum is over all classical trajectories that go from q' to q on PSS and correspond to 
one Poincare mapping. All trajectories should leave the PSS with the positive, say, normal 
component of the momentum. The reduced action S(q,q';E) is still measured in the units of 
h. The Maslov index v is determined by the number and type of caustics encountered by the 
trajectory pi). For example, a regular caustic of dimension (N — 1) increases v by 1; a hard 



wall increases v by 2. We discuss the origin of Maslov phase in Sec. 1.6. In a billiard with its 



boundary as a PSSQ v = 2 and there is only one trajectory q' i — ► q. Then Eq. (1.23) follows. 
A surface of section wavefunction is determined from the integral equation 

(1.31) 



V(g)= / dq'T(q,q')i>(q') 
JPSS 

that we will call Bogomolny l 's equation and the quantization condition is 



det(l — T) = 0. 



(1.32) 



The original wavefunction can be found by propagating ip(q), 



tf(r) = / dqG(r,q;E)i/>(q) 
JPSS 

where 



G(r,q;E) 



(2wi) 



(N-l)/ 



l)/2 E 



cl. tr. 



v(r) 



det 



d 2 S(r,q;E) 



dr±dq 



1/2 



JS(r,q;E)-i%v 



(1.33) 



(1.34) 



In the billiard problems it is sometimes convenient to choose PSS infinitesimally close to the hard wall, not 
on the hard wall. 
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The trajectories start at point q on PSS with the positive normal component of the momentum 
and end at point r of the interior; v(r) is the modulus of the velocity at point r; r± is the 
direction perpendicular to the trajectory at point r. In the case of a billiard one can show that 



Eqs. (1.33) and (1.34) agree with Eqs. (1.10) and (1.21) up to normalization if the asymptotic 
form of HankcPs function is used. 



1.5 Properties of the T-operator 

In this section we review some of the properties of the T-operator p9[ |. As was noticed by 
Bogomolny, T-operator is similar to the time-dependent semiclassical Green's function [cf. Eq. 



(1.25)]. It is therefore not surprising that they share some of the properties. To simplify 
notation we assume the system to be two-dimensional and omit the Maslov indices. 

First, note that the T-operator vanishes for short trajectories, i.e. T(q,q') — > as q — > q' . 



Indeed, in this case S(q, q') ~ p \q — q'\, so the prefactor in Eq. ( 1.30| ) vanishes 



As was mentioned before, the T-operator, unlike the exact kernel K(q,q'), is unitary in 
semiclassical approximation. To see this we evaluate the integral 



d 2 S(q, qi )d 2 S(q,q 2 y 1/2 



e <[S(9,9a)-S(9,9i)] (1.35) 



dqdqi dqdq2 

by the S$. The stationary phase condition dS(q, qi)/dq = 0S(q, q 2 )/dq requires that the orbits 
<7i i — ► q and 172 1 — ► q have the same final momentum in the PSS. For a deterministic map it 
means they are the same orbit, i.e. q\ — q 2 . Expand S(q, q 2 ) — S(q, qi) + (q 2 — Qi)9S(q, q\)/dq\. 
Then the integral becomes 

(TtT) (q u q 2 ) = ±- j " dq e^fe-*) s S (q 2 - qi ). (1.36) 

The sum over classical orbits disappeared because the integration over dpi takes care of all 
the orbits leaving q\ . The 8- function is not ideal since \p\ \ is limited by the fixed total energy 
but it is a good approximation semiclassically. The unitarity of the T-operator leads to the 
resurgence of the spectral determinant J|3| . 

Another useful property, the nth power of the T-operator has the form of the T-operator 



for n Poincare mappings, that is T n (q,q'; E) is given by Eq. (1.30) but the orbits now cross 
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the PSS n times. Consider, for example, 



The 5$ condition dS{q,q")/dq" + dS{q" \q')/dq" = ensures that the orbits g' i — ► q" and 
g" i — > q have the same momentum at q" . In other words, it selects the classical orbits g' i — > g" 
i — ► g that cross the PSS twice. To complete the proof we need to show that 

d 2 S(q,q") d 2 S(q",q') 2 

a 2 g2(g,.?') <9g"<9g' 1 j 

9g" 2 

where $2(5, <?') = ^(Sj + <S(g", g'). For this we differentiate the condition by q and by 
q' . The composition property is important for the derivation of the trace formula using the 



T-operator (see Sec. 6T) and for the semiclassical Fredholm theory |33|, |3J, |3C[| . 

The T-operator in momentum representation is semiclassically a finite matrix if the PSS 
has a finite length C. Suppose the PSS is a closed line. We define a complete set of the 
momentum eigenfunctions 

Mq) - 4=e^ (1.39) 
where p = (2ir / C) x integer. The T-operator matrix is 

T PP . = I / dqdq'T(q,q'y(p'l'-™). (1.40) 

The 5$ condition requires that p and p' be equal to the projection of the classical momentum 



at points q and q' , respectively. But this projection is limited by y2mE by the absolute value. 



Hence the size of the matrix is no greater than C\f2mE / V or 2£/A where A is the de Broglie 
wavelength. 



1.6 Maslov phase and the change of variables 

A neighborhood of a caustic in a multidimensional system, like a turning point in one dimen- 
sion, is the region where the semiclassical approximation breaks down. It turns out that the 
T-operator and the wavefunctions acquire additional phases, called Maslov's phases, in the 
regions of validity pi] . The Maslov phases will change, in general, the overall probability 
density and the quantization conditions and thus it is important to understand their origin 
and be able to properly account for them. 
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The central point in Bogomolny's derivation of the T-operator is the use of the semiclassical 
approximation to the Green function in the energy representation |40j 



G(r,r';E) - - — {N _ 1)/2 



cl. tr. 



det 



d 2 S{r,r';Ey 1/2 



v(r)v(r') 

JS{r,r':E)-i%i> 



dr±dr' ± 



It satisfies (E — H)G(r, r'\ E) — 5(r — r') and is not valid for short trajectories, at least in two 
dimensions. rj_, r'j_ arc the local coordinates perpendicular to the velocity at r, r', respectively. 
v = \dH/dp\ is the modulus of the velocity (actions and momenta are defined in the units of 
h). The Maslov index v in the T-operator ( 1.30| ) is inherited from this Green's function. 



Let us follow a particular orbit starting from r' . At some points along the orbit one or 
more eigenvalues of the matrix d 2 S/drj_dr' ± may become infinite and change sign pT[ . That 
is what happens when the orbit touches a caustic (see below). If m eigenvalues change sign 



(caustic of order m) the determinant in Eq. (1.41) changes by (— l) m . If we drop the modulus 



and analytically continue the determinant in the complex space "around" the singularity |4q| 
as in the one-dimensional case |f9| we will find that the Green function acquired the phase 
— m7r/2. One has to assume the function is exponentially small in the classically forbidden 
region beyond the caustic. (This method does not work near a hard wall.) The Maslov index 
of an orbit is the sum of the contributions from all caustics the orbit encounters. Each caustic 
increases v by m, each hard wall increases v by 2. In particular, if r' — q' and r = q lie on the 
surface of section and the orbit makes one Poincare mapping, this index enters the T-operator. 

So far we defined the Maslov index in terms of singularities of the matrix d 2 S/ dr±dr' ± . Now 
we give it the geometrical interpretation j2^, [H| . Let us surround the orbit of the previous 
paragraph (we call it the central orbit) with a sufficiently narrow tube of trajectories with 
the same energy E all of which originate from point r'. They form an (TV — l)-parameter 
family that can be parametrized by the vector dS/dr^, which depends on r. For a given orbit 
in the family we will choose r to be on the caustic. It is possible that all the trajectories 
within the tube (remember, the tube is narrow) touch an (N — m)-dimensional surface which 
is independent of the shape of the tube. This surface is called a caustic of order m. Suppose 
the central orbit touches the caustic surface at point r. The caustic has m — 1 dimensions less 
than the number of parameters in the tube. This accounts for m — 1 infinite eigenvalues of 
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t(r) Ax 



caustic 



- (Ax) 




Figure 1.2: Two orbits that touch the caustic and lie almost in the same plane. 



d 2 S/dr±dr' ± . (When different trajectories touch the caustic at the same point, the derivative 
of the parametrization vector becomes infinite.) Without the loss of generality we can always 
assume that the caustic is locally flat by making an appropriate coordinate transformation. 
The orbits are parabolic near the caustic. Consider an orbit that touches the caustic a small 



distance Ax from the central orbit and lies locally in almost the same plane (Fig. 1.5). (To 
find this orbit note that the plane of the parabola of the central orbit is fixed by two unit 
vectors: one is the tangent vector t(r) at the point of touch, another, n(r), is orthogonal to 
the caustic at the point of touch. If the neighboring orbit touches the caustic at r + t(r)Ax, its 
tangent at this point is almost parallel to t(r) by continuity. Its second vector n[r + t(r)Ax] 
belongs to the m-dimensional complement to the caustic, so its direction can be fixed by m— 1 
free parameters mentioned above.) Now, the difference in parameters for these two orbits is of 
order Ax, but the distance between them at point r in the direction of n(r) is of order (Ax) 2 . 
This brings another singular eigenvalue of d 2 S/dr±dr' ± . Thus, this matrix has m singular 
eigenvalues in total. 

A special case is the caustic associated with the original point r' . Creagh et al. ]2(j] argue 
that each negative eigenvalue of —d 2 S/dr±dr' ± for positive time t = 0+ decreases the Maslov 
index by —1. (Here one has to be careful to preserve the orientation of the local basis (ey , e^_), 
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i = l,...,N — 1, along the orbit, where ey points in the direction of propagation.) This 
apparently contradicts the statement in Ref. jll| that v = for short trajectories since there 
are no caustics in between. Notice, however, that for the kinetic-plus-potential systems all 
eigenvalues of —d 2 S/dr±dr' ± are positive (^6). In general, as the example in the next section 
shows, this additional Maslov index is important for the quantization conditions. 

It is sometimes convenient to evaluate the classical action in the canonical coordinates, other 
than the original ones. For example, in the next section the action-angle variables are used. 
It would be desirable to be able to construct the T-operator and solve for the wavefunction 



directly in the new coordinates. Let us analyze how the Green function (1.41) transforms 
under canonical transformation (p, q) i — > (P, Q). We rewrite Eq. ( p. .41 ) in an equivalent form 


y/A(q, piS{q , q :, E ) 



G(q,q';E) = £ G 1 (q,q';E) = £ A___ £r 

cl. tr. 1 y^Ki) 



cl. tr. 



*M:*> = -(0) 



1-N 



dct 



(d 2 S d 2 S d 2 S d 2 S 



(1.42) 



(1.43) 



\dE 2 dqdq' dEdqdEdq' , 
Note that we removed the modulus from under the square root together with the Maslov 
index. For a given trajectory we can define the coordinate transformations Q — Q(q, q') 
and Q' = Q'{q,q') which can be multivalued. The probability density for a semiclassical 
wavefunction is a sum of probabilities for individual orbits, therefore the transformed Green's 
function 

G^\Q 1 Q')= ]T G x {q,q')\ 



cl. tr. 



det 



dQ 



det 



dq' 



dQ' 



where q and q' are the functions of Q and Q' . But A(q,q') 
det (dQ' /dq') . Finally, the Green function 



(1.44) 

A(Q,Q') det (dQ/dq) 



GW(Q,Q';E)= ]T JA(Q,Q';E) sgn ( det^det^ 



cl. tr. 



^iS {Q) (Q,Q';E) 



dq dq' J i (2m) {N - 1)/2 



(1.45) 



where S^(Q,Q';E) — S(q,q']E). Equation ( 1.45 ) has the same form as Eq. ( 1.42 ) except 
for the signature of the Jacobians. Let us fix q' and follow q along the trajectory. If Q — Q(q) 
is a multivalued function, some eigenvalues of dQ/dq become singular and change sign as 
the orbit goes from one sheet of Q(q) to another. This happens if a caustic in q-space is 
completely or partially removed in Q-space by making Q(q) multivalued. If, for instance, 
the caustic is removed completely, the number of singular eigenvalues is equal to the order 
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of the caustic. Suppose k eigenvalues become singular. Then the matrix inside A(Q,Q';E) 
has k singular eigenvalues less then the matrix inside A(q,q' ; E). Thus the total number of 
singular eigenvalues, and, consequently, the Maslov index, is invariant under the coordinate 
transformation.^ The important conclusion from this discussion is that the semiclassical Green 
function or the T-operator can be evaluated with Eqs. (1.41) or ( 1.3C| ) in any set of canonical 
coordinates, apart from the Maslov index which is determined by the topology of the orbits in 
the physical space. This result is quite different from the classical assumption that all canonical 
coordinates are equivalent. 



1.7 Two-dimensional separable system 

As an illustration of the T-operator technique we apply it to a two-dimensional system which 
is separable in the Cartesian coordinates. We will solve the problem in the action-angle (A A) 
variables. They change the topology of the orbits and thus bring certain complications. It will 
be helpful to understand them before moving on to the perturbation theory of the following 
chapters. The results are, of course, well known from simpler methods. The quantization of 
the rotationally invariant integrable systems in two and three dimensions using Bogomolny's 
T-operator was done by Goodings and Whelan p6|| . 

For a separable system we can define the action-angle variables (I X ,9 X ) and (I y ,9 y ) in x 
and y directions, respectively. We choose the x axis (y = 0) as the PSS which corresponds to 
Q y = const. One Poincare mapping x' i — ► x is described as 0' x i — ► 9 X , 9 y i — ► 9 y + 2ir in the AA 



variables (Fig. 1.3). Note that by definition 9 X > 9' x ; the actions ij = (2n)^ 1 <j> pidxi (i = x, y) 
are integrals of motion; the functions x{9 x ) and y{9 y ) have period 2n. We made a transforma- 
tion to the AA variables because in these coordinates the T-operator has a particularly simple 
form: 



T{9,9'-E) 



1 

2ttI 



d 2 s 



1/2 



e iS(Ae;E)-ifu(e,e';E) ( 14 g) 



d9d9' 

where we dropped the subscript "x" in 9 X and, apart from the Maslov index, T(9, 9 '; E) depends 
only on the difference A9 = 9 — 9'. We assumed there is only one trajectory connecting 9' 

2 If there is a new caustic of order m in Q-space that did not exist in q-space then m eigenvalues of dQ/dq 
go through zero and change sign but in the opposite direction than the corresponding m singular eigenvalues 
of A(Q, Q'\ E). Thus, again, the Maslov index does not change. 
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Figure 1.3: Poincare map from 9 X = 0' to 6 X = 6 on the surface of section 9 y — const. 



and 6. If there are more than one trajectory with the same energy, they will have different 
I's, and the PSS wavefunctions ip (0) cx e llmd can be treated independently in the semiclassical 
approximation. 
The action 

S(A6; E) = A6I X + 2i:I y (1.47) 
where the 7's should be expressed in terms of A8 and E. One relation is the Hamiltonian 

H(I x ,I y ) = E. (1.48) 

Suppose it can be solved for I y , 

Iy=9E{I x ). (1.49) 

Define the frequency ratio 

u x (dH/dI x ) Iy 



where I y was substituted from Eq. ( 1.49| ). Then the second relation is 



a E {I x ) = (1.51) 

Z7T 



Eqs. ( 1.49| ) and (1.51) allow one to write the action S(A9; E)^ Note a useful formula j76| 
a = -g' E (I x ) (1.52) 

3 The T-operator is not well defined for a two-dimensional harmonic oscillator where A8 = const. 
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which can be derived by differentiating Eqs. ( 1.4S| ) and ( |l.49| ). 

The Maslov index v(9, 9'; E) counts the number of caustics encountered by the orbit in the 
Cartesian coordinates. Suppose there are two turning points in both x- and y-directions that 
produce the first order caustics. The Maslov index 



v(e,6') = v + v x (6,9') + u y) v y = 2. 



(1.53) 



Here v x , v y are the number of times the orbit goes through the turning points in x- and 
y-directions, respectively, and Uq = —®(d 2 S/drxdT , _ L |t_>o) is related to the starting point 
caustic (see the previous section), O is the Heaviside function. Assume for the definitiveness 
that the turning points in x are located at #mod27r = 0,7r (they must be separated by 7r 
because of the time-reversal symmetry). Suppose Ttn < A9 < 7r(n + 1), for some n = 0, 1, . . . . 
Then, depending on the end point 9, the orbit encounters n or n + 1 turning points in x. It is 
easy to find 

n + l, < 9 < A9 - nn 

(1.54) 

n, A9 — im < 9 < tt 

continued in 9 with period ir. The additional Maslov index vq can be conveniently expressed 
as © 

= -6 . (1.55) 



u x (9,A9) 



To solve Bogomolny's equation (1.31) we make an ansatz 



- p ue+if{9) 



(1.56) 



where / = const to be determined and f{9) is a step-function. With T given by Eq. §LA(\ ) 
Bogomolny's equation becomes 



Jf(8) 



d{A9) 



\S" (A9)\ 



2m 



1/2 



1 i[s(Ae)^/Ae+/(e-A0)-f v(e,AS)\ 



(1.57) 



where the integration variable has been changed from 9' to A9 and e lie was canceled on both 
sides. The stationary phase condition is 



I = S'{Ae st ) = I x (A0 st) E). 



(1.58) 



The second equality follows from Eq. ( 1.47 ) taking into account that ( j = g' E {I x ) ( ) 



= — ^ (§as) E [ see Eqs. ( 1.51 ) and ( 1.52Q ]. Note that / and v do not change the stationary 
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point. We expand S (A9) near the stationary point and integrate. Equation ( 1.57 ) becomes 

|S"(A(? st )p 1/2 



S"(A9 st ) 



o i[27r/ y (A0 st ,S)+/(0-A0 st )-fi/(0,A0 st )] 



(1.59) 



Applying another chain of equalities S" {A9) — (-§^0) E = 57 (^) E = — 
the pre-exponential factor cancels the Maslov index vq. Let us require that 



f(9) - f(9 - A9 st ) + - Vx {9, A9 st ) = 



„ 1 „ wc find that 



(1.60) 



for all 9. Then Eq. (1.59) is solved if 



I y (A9 st ,E) =n y + 



(1.61) 



for some n„ 



0,1,. 



This is, of course, a well known Einstein-Brillouin-Keller (EBK) 



quantization of an action variable. One can check that the step-wise function (Fig. |l.4|) 

f(9) =0if 0<9<tt; f(B + «)-f(6)=~ (1.62) 



satisfies Eq. (1.60). The second quantization condition 

I x {A6 st ,E) =n x + ~, n x = 0,1,. 



(1.63) 



comes from the 27r-periodicity of ip (9) Eq. (1.56). (This requirement ensures that ip (9) is 
uniquely defined on the cross-section of an invariant torus < 9 < 2n, 9 y — const.) 
Thus, after all this trouble, we arrive to the familiar EBK expression for the energy 



E, 



n. ,n y 



H I 



2 ,Iy 



with the PSS wavefunction 



4,(0) 



(1.64) 



(1.65) 



This function is discontinuous at 9 — nn. This may seem to be unphysical, but remember that 
the semiclassical approximation breaks down near the turning points. Therefore Eq. ( 1.65| ) 
should not be used near 9 = im. The exact wavefunction would smoothly join the discontinuity. 
If there were hard walls instead of caustics the actions would be integral and f(9) would have 
the steps of size tt. 

The wavefunction in ^-representation can be recovered by the standard procedure ip (x) = 



ip (9) yj\d9/dx\. Note that 89 jdx = lu x /v x , where v x is the velocity, and I x 9 = J p x dx. Suppose 
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Figure 1.4: Step-wise phase f(8) from Eq. ( L5q ) 



x{6 = 0) is the turning point on the left. Then ip (x) has two branches: 

^ ( x ) oc _1 e */" ( o)lP-l«te' ) O<0<tt, 



'f'.r 



^ 2 ( a! )«-i=e- iJ '-w l, '- |,to ' +<f , -7r<0<O, (1.66) 



which makes 



V' (a;) oc — == cos 



the one-dimensional bound state wavcfunction |4 



/ |p x |da;-i-- 

Jx(O) 4 



(1.67) 



Finally, we note that the algorithm we used to solve the semiclassical problem in the 
AA coordinates (or other coordinates that change the topology of the trajectories) may be not 
unique. Both in this and in the previous section, where we considered the change of coordinates 
in the T-operator, we assumed that the wavefunctions do not acquire Maslov's index on caustics 
— we always took the absolute value of the Jacobian: ip (g) = ip (Q) \J\dQ/dq\ . Alternatively, 
we could drop the modulus or add the Maslov phase. Then the Maslov phase of the T-operator 
would be determined by the caustics in the current coordinates (not physical space) and the 
wavefunction ip would not have the unnatural discontinuities as in the above example. However, 
the straightforward quantization conditions would be incorrect, for instance, the actions would 
be integral instead of half-integral. Thus the special formulation of the quantization conditions 
would be necessary, again giving significance to the physical coordinates. Whichever method 
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is used, the physically important quantities, the wavefunction tp (x) and the quantized energy, 
will be the same. 

1.8 Conclusions 

Bogomolny's T-operator is a powerful tool that combines the scmiclassical approximation and 
Poincare's surface of section. Bogomolny's equation determines the surface of section wave- 
function and the energy levels. The full wavefunction can be reconstructed from the surface 
of section wavefunction with the help of a semiclassical propagator. In billiards the method 
is related to the boundary integral method. The properties of T-operator include unitar- 
ity, the composition property, and the finite size in the momentum representation, all within 
the stationary phase approximation. The T-operator propagates the semiclassical wavepacket 
consistently with the classical surface of section map. The regions where the semiclassical ap- 
proximation breaks down, like caustics and walls, are responsible for the additional (Maslov) 
phases in the T-operator and wavefunctions. One should be careful, when making coordinate 
transformations with singular points, not to change the Maslov phase. For the integrable 
systems the T-operator method is consistent with the EBK quantization. 
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Chapter 2 



Perturbation theory 



In this chapter we present a systematic derivation of the perturbation theory that provides a 
semiclassical description for almost integrable systems |}2) . In a theory of this type there are 
two competing quantities: the large action (compared to the Planck constant) and the small 
perturbation from integrability. Their interplay determines to what degree the structure of 
classical phase space is reflected in the quantum results. 

It is a well-known fact of the classical theory that the topology of invariant tori near 
the periodic orbits changes under perturbation. One needs a resonant classical perturbation 
theory in order to describe it. On the other hand, the standard EBK quantization procedure 
relies on the classical invariant tori, which means that a resonant semiclassical theory might be 
necessary for some classes of states. Here is a simple estimate. For a periodic orbit of action 
I the perturbation of size e changes the topology of the tori within the layer SI ~ \fel. If 
51 > H, several quantum levels become mixed, i.e. the resonant theory is required. 

For a billiard of linear size L with the perturbed boundary the condition becomes kL^/e > 1 
where k — 2ir/\ is the wavenumber. This means that, unless the perturbation SL -C A 2 /L, the 
billiard will have strongly perturbed states .0 Note that even the distortion of the boundary 
shorter than the wavelength may be strong enough to mix several energy levels. 

In the following sections we derive the analytic expressions for the wavefunctions and energy 
levels of the states associated with the classical resonances. Although our theory is effectively 

a billiard the level spacing h 2 /mL 2 should be compared to the shift of a given level h 2 k 2 SL/mL, which 
results in the same criterion. (Noted by M. Sieber.) 
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a resonant perturbation theory, the diagonalization of unperturbed states is not explicit. The 
Bogomolny equation allows to express the results in a simple, easy to visualize form. 

We use the circle billiard with a perturbed boundary as an example. Several papers pub- 
lished in recent years plj , |32| , p0| , p5| discuss the localization and diffusion in angular momen- 
tum space of this system. We too find the angular momentum localized, with the degree of 
localization depending on the smoothness of the perturbation. 

There are several other methods that deal with almost integrable systems that are some- 
times similar, sometimes complementary to the following theory. For instance, the quantization 
of Birkhoff-Gustavson normal form |L5], |3^, ^38| is useful for the perturbed harmonic oscillators 
— the case to which our theory does not directly apply. The Born-Oppenheimer approxima- 
tion (Ch. ^|) is similar to the perturbation theory in some cases. The perturbed Berry- Tabor 
formula (Ch. ^) expresses the contributions of the periodic orbits to the density of states 
consistently with the perturbation theory. 



2.1 Perturbed integrable systems 

An ./V-dimensional Hamiltonian system is called integrable if it has N independent integrals 
of motion in involution fr4| . Its motion in the phase space is confined to an .ZV-dimensional 
manifold that has the topology of a torus. It is called an invariant torus. The integrals of 
motion can be chosen to be the action variables Ij, i = 1, . . . , N. Then the angle variables 9i 
(identified with 8i + 2tt) provide the natural coordinates on the torus. The Hamiltonian of the 
system depends only on the actions, H (I) = E. The equations of motion Oi = dH /dli = uji(I) 
can be easily integrated. 

In the case N = 2 we can choose the Poincare surface of section (PSS) as, say, 82 = (Fig. 
|l.3| ). The dynamics of the system induces the PSS map (I',0') 1 — ► {1,6) ( we dropped the 
subscript "1") simply as / = /' , 9 = 9' + 2tiuoi/lo2- If all points (1,9) belonging to a certain 
trajectory are plotted I vs 6, < 6 < 2tt, they all will lie on a horizontal line I = const. The 
line will either be covered densely if u>i/u>2 is irrational or have a finite number of points if it is 
rational. We may say that the line / = const is an intersection of an invariant torus with the 



PSS. The reduced action S(6, 9') = S (6-6') = (6-6')I 1 +2irI 2 [cf. Eq. (1.47)] is a generating 



function for this map: I = dS/d9, I 1 = -dS/86 1 (Sec. 1.3). Note that S(9,9') depends only 



20 



on the difference 9 — 9'. The action variables must be regarded as functions of 9 — 9' and E 
(see Sec. [L7| for details). 

The Hamiltonian of a perturbed integrable system 

H = H {h,h) + eH 2 (h, la, 01, 02 ) + e 2 H 4 + ■■ • (2.1) 

differs from an integrable Hamiltonian by the terms proportional to a small parameter e. 
Likewise, the action of such system 

S(9, 9') = S (9 - 9') + eS 2 (9, 9') + e 2 S± + --- (2.2) 

is expanded in powers of e. (We reserve the odd subscripts for half-integral powers of e that will 
appear in other expansions.) The classical perturbation theory [p2| , in principle, allows one 
to find the corrections to the action if the Hamiltonian is known, for example, 5*2 = — J H 2 dt. 
Since only short orbits are involved, the problem of small divisors does not arise. In the case 
of a billiard with a perturbed boundary, S(9,9') can be easily deduced directly. 

Of course, when the perturbation is on, Fs are no longer integrals of motion. Consequently, 
the orbits on the PSS are no longer confined to the lines / = const. One can imagine several 
possibilities. The lines could be slightly deformed, they could change their topology, or disap- 
pear completely (i.e. the orbit would densely cover a two-dimensional area in the PSS). As a 
matter of fact, all three cases may be present within one system, depending on the region in 
the phase space. 

According to the Kolmogorov-Arnold-Moser (KAM) theorem, for small enough e and 
smooth enough perturbation most of the tori continue to exist and remain close to the unper- 
turbed tori. The measure of the destroyed or considerably modified tori tends to zero as e — * 0. 
The destroyed tori are located near the rational, or resonant, original tori. These are the tori 
that have a rational winding number u)\/uj2, i.e. they support a periodic orbit, or resonance. 
Consider a (p, q) resonance, i.e. the winding number is p/q, an irreducible fraction. Then all 
points on the line / = const are the fixed points oiT q , where T is one Poincare mapping. The 
Poincare-Birkhoff fixed point theorem j^J states that under the perturbation only an even 
number, a multiple of 2q, of fixed points will remain. They will alternate between stable and 
unstable. 



As Fig. 2.1 illustrates, the tori form resonant islands around the fixed points. The size 



21 



of the islands is proportional to y/e. Thus most of the tori near a rational torus change their 
topology under the perturbation and form a new system of tori. The new tori, in turn, form 
secondary resonance chains (as in the example of (2, 5) resonance in the figure) although they 
are exponentially smaller than the primary resonances [ |52|| . There is a chaotic region near 
the separatrix formed by numerous intersections of stable and unstable manifolds. It too is 
exponentially small p0[. The size of the resonant islands reduces rapidly with q (see Sec. 



2.3.2), so the total measure of the strongly modified phase space is finite and proportional to 



y/e. The rest of the tori are slightly perturbed on the scale of e. Figure 2.1 shows the torus with 
the "most irrational" golden mean (GM) winding number. It is supposedly the last torus to 
be destroyed with the increase of the perturbation, although this is not a subject of this work. 
The figure also illustrates how the rational numbers with large q tend to become "closer" to 
irrational numbers. If one disregards the small island chains in (7, 16) or (15, 32) resonances 
the overall curve is similar to a perturbed irrational torus. Its oscillation is of order of e. This 
effect is relevant to the semiclassical theory since the quantization misses those details of the 
phase space that have a typical area less than 2ttH (see Sees. 2.3.1 and 2.3.2] ). 

If the perturbation is not smooth enough (so that the KAM theorem does not apply) the 



invariant tori may not exist at all [Fig. 2.2 (a)]. In this case the semiclassical perturbation the- 
ory needs diffraction corrections, although the simplest version still reflects the main features 
of the states for small e. 



2.2 Formulation of the theory 
2.2.1 Perturbed circular billiard 

The perturbation theory we are about to present is based on the T-operator and requires 
a careful choice of coordinates and surface of section. From the theoretical point of view 
the natural choice would be the action-angle (AA) variables where the unperturbed action 
depends only on the difference of the angles. In practice working in the AA variables is 
not always physically transparent and may require a complicated transformation from the 
original coordinates and back. On the other hand, billiards with a perturbed boundary have 
an advantage that the reduced action is just a geometrical quantity proportional to the length 
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Figure 2.1: The surface of section (boundary phase space) of a perturbed circular billiard, l/k 
is an angular momentum normalized to its maximal value and 9 is an angular coordinate along 
the boundary. The numbers on the right denote the orbits near the (j>, q) resonances perturbed 
by order ^/e. GM is an orbit with the golden mean winding number perturbed by order e. 
The secondary resonances and the chaotic regions near the separatrices, that are exponentially 
small when e — > 0, can also be observed. Here the perturbation is the "smoothed stadium" 
with the parameters e = 0.1, r? = 0.19 defined in Sec. 2.2.1. 
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Figure 2.2: The surface of section for (a) "short stadium" and (b) "smoothed stadium" of Sec. 
2.2.1 with e — 0.01, rj = 0.15. l/k as in Fig. 2.1. The points O, +, and * denote the intersection 
of the invariant tori with the surface of section, h nv (6), as given by the perturbation theory (see 
Sec. 2.3.1), for the non-resonant, separatrix, and resonant orbits, respectively. In (a), where 
KAM fails, orbits started at symbols were iterated forward and backward 15 times (appear as 
dots •). Only short time structure is regular. In (b), where KAM applies, three orbits, each 
iterated 1000 times, coalesce into solid lines, (c) Husimi plots for the numerical states of Fig. 
2.4, as well as the scar state V'32, for the short stadium. The square has area h. 
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of the orbit. A perturbed circular billiard whose boundary in the polar coordinates is r{9) = 
Rq + eAR(9) is a particularly convenient system since the angle 9 together with the (modulus 
of) angular momentum I are, in fact, the AA variables (when e = 0, of course). (We can 
always make the angular average of AR vanish and assume AR ~ Rq.) We will formulate the 
perturbation theory for this system, which allows the simple mathematical description and 
direct physical interpretation. There should be no principal difficulty to generalize the theory 
or tailor it to other systems, as we do in the subsequent chapters. We would like to point out 
some limitations that will become more apparent as we proceed. First, the theory is the most 
useful in two-dimensional systems, where the PSS has a one-dimensional space component. 
Second, the unperturbed action So should depend only on the difference of coordinates in the 
PSS which may force one to use the AA variables. Finally, the theory is not suitable for the 
perturbed two-dimensional harmonic oscillator since So is not well defined (see the footnote 
on p. 0). 

From now on we assume Rq = 1. A wavenumber k also becomes dimensionless. (To restore 
the proper units one should substitute k by kRo.) A perturbation of the form AR(9) = 
| sin 9 1 — 2/n models the short Bunimovich stadium that recently received a lot of attention 
|(J ||. The straight segments in this stadium have length 2e. Our perturbed circle 



approximates it to the order of e and has a discontinuous derivative (Fig. 2.3) [we neglected 
the straight segments when deriving AR(9), which thus describes the outer boundary of two 
circles with their centers being 2e apart]. Although the invariant tori do not exist, our theory 
still works if fee 3 / 2 <C 1. We also consider a "smoothed stadium" AR(9) = \J sin 2 9 + rj 2 — C v , 
where C v makes the angular average vanish. If e is sufficiently small compared to rj the 



invariant tori exist [Fig. 2.2 (b)]. We choose rj ~ y/e in our numerical examples. 



We associate the PSS with the boundary of the billiard. We parametrize it by the angle 



9 with the conjugate momentum I. Our goal is to construct the T-operator ( 1.23 ) and solve 



Bogomolny's equation ( 1.22 ) perturbatively. For convenience we assume H = 1 and the mass 
m = 1/2. Then the energy E = k 2 . The action 

S(9, 9') = kL(0, 9') = k [L Q {9 - 9') + eL 2 {8, 9') + ■ ■ •] (2.3) 

where L(9,9 r ) is the chord length between points 9 and 9' on the boundary. For a perfect 
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Figure 2.3: Short Bunimovich stadium R(8) = R + eAR(6) with e = 0.3. The straight 
segments have length 2e. The dashed line denotes the circle of radius Rq = 1. The difference 
between the actual stadium and the approximation for AR(6) that we use is of order e 2 . The 
classical orbits that pass through the shaded region are affected by the I = resonance. 



26 



circle wc find 



and the correction 



L (B - 9') 



Ln 



[AR(9) + AR(8')} 



(2.4) 



(2.5) 



We always assume k 3> 1 and e <C 1 as general requirements for the semiclassical and pertur- 
bation theories, respectively. 



An example of the PSS is shown on Fig. 2.1. Note that in a perfect circle of radius 1 the 



angular momentum is bounded, |/| < k. The resonances correspond to the periodic trajectories 
in a circle. For instance, the (1, 2) resonance is related to the I = orbit passing through the 
center of the circle. The (1,3) resonance is a triangle shaped periodic orbit, etc. Our theory 
is designed to find the quantum states located primarily within a given resonance. 



2.2.2 Simplest case: (1,2) resonance, fce < 1 

The theory simplifies significantly for the states near (1,2) resonance, partly due to the sym- 
metry of the perturbation L2(0, 9'). We also assume the perturbation to be sufficiently weak, 
ke <C 1, a requirement that will be later relaxed. It was explained in the chapter's introduction 
that as long as ky/e > 1 the resonant perturbation theory is needed. The relevant classical 
orbits have the angular momentum I ~ fc^/e and pass through the shaded region in Fig. [2.3] . 

In a perfect circle the resonant trajectory maps point 9' to 9 = 9' + tt. We therefore expand 
the T-operator in 59' = 9' — 9 + tt, as well as in ke. In this case Eq. ( 1.23j ) becomes 



T(6, 8') ~ -J— (1 + • • .)e i ( 2fc -i ae ' 2 +-) [1 + ikeV(8) + • • •] (2.6) 

V 47TI 

where V{6) — L 2 {9, 9 — tt) = AR(8) + AR(9 — tt). We also expand the wavefunction 



^(0') = i]){6 - tt) + SO'tp' 



q/2 



-^"{9 - tt) 



(2.7) 



that solves Bogomolny's equation ip = Tip. We evaluate the integral J T(Q, 9') %p(d')dd' in 
the stationary phase approximation (5$) thus justifying the smallness of 59' . The equation is 
reduced to 



ip(6 -tt)- -ib"(6 -tt)+ ikeV(9U(9 - tt) 
k 



(2.8) 
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Assume we can expand ip(6—ir) ~ e~^ 2k+ ^ ) i/)(0) — ikeE m ip(6—n), making i/j(9-tt) = e~ iu ip(9) 
where 

cj(fc) ~2fc + | + /ce£: m . (2.9) 

The constant £" m is to be determined. We find that ip(9) has to satisfy the ordinary differential 
equation 

4>"(e) + k 2 e [E m - V(9)} 1>{ff) = 0, (2.10) 

which is a one-dimensional Schrodinger equation. We will refer to V{9) as a "potential" and 
to E m as a "(surface of section) energy" when it does not cause a confusion. Note that the 
strength of the potential is scaled as (fci/e) 2 . 

The eigenstates ip m (9) of Eq. ( 2.1 0| ) with the eigenenergy E m solve Bogomolny's equation, 
m is a quantum number that labels the states. Since V{9) has period 7r we can always find 
the 27r-periodic solutions such that ip m (9) = ±ip m (9 — n). The second quantization condition 
is then 

uj(k) = 2im + oj m (2-11) 

where w m — or 7r (depending on m) and n is integer. This provides an equation for k and 
the total energy E nm 

, 2-7T (n — + u> m 

knm — V Aim — ~ ^ ■ (2-12) 

Note that E m and ip m weakly depend on n via the parameter k 2 e. 

Eq. ( 2.1 0| ) can be solved by standard methods either numerically or analytically. If k^/e 3> 1 
we can approximate its solution in WKB: 

ex -L=e* k ^e) ( 2.i3) 

where 

= ± J d6y/E m -V(0). (2.14) 
If £J m > V(9) for all 0, we may quantize E m by the condition 

ky/e [ d9^Em - V{6) = u m + 2%m. (2.15) 
Jo 

We call this the "rotational" case since all 9 are in the "classically allowed region." As will 
be clarified later, these states quantize the classical orbits in the surface of section that lie 
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near the resonance chain though outside of the separatrix (Fig. 2.1). In the KAM language, 
the topology of the tori did not change as a result of the perturbation, but they are strongly 
distorted by the nearby resonance. If we compare a picture of a resonance chain with the 
phase portrait of a pendulum, these tori would be analogous to the rotational trajectory of a 
pendulum. 

Likewise, if E m < max V(9) the motion will be librational. These are the tori that lie 
inside the separatrix and surround the stable periodic orbit. Their topology is changed by the 
perturbation. The potential V(9) has at least two wells between and 2ir. If the tunneling 
between the wells can be neglected the degenerate levels will be given by the Bohr-Sommcrfcld 
conditions p9[ 

= 7r (^m + — \ . (2-16) 

Here the limits of the integration are the classical turning points. The wavefunction inside the 
well is 

ij m (e) = lE m -V(6)}- 1/4 sm( f d9'^E m -V(6') + j\ . (2.17) 



Figures 2.2 (c) and 2.4 show the examples of the rotational and librational states. The two- 



dimensional wavefunctions for the librational states with m = 2, 3 are shown in Figs. 2.5 
respectively. 

2.2.3 General case: (p, q) resonance, ke'l 2 <C 1 

In this section we quantize the orbits near an arbitrary (p, q) resonance. We also lift the re- 
striction ke <C 1. Instead, we construct the perturbation series ip(9) — exp {ik [fo{9) + bfi(9)+ 
b 2 f2{9) + ■■■ + b M f M (9)] } , where b = and f[ ~ 1. We call M the order of the perturba- 
tion theory. The large and small parameters k and b should satisfy the conditions ^> 1 



and kb M+1 <C 1. Equation ( 2.13 ), for example, is the case M = 1. The first condition ensures 
the validity of this WKB-type ansatz by making the phase fast varying. In other words, the 
quantum system will reflect the classical structure to order M. The second condition allows to 
neglect the terms of order M + 1. In practice, the procedure becomes very tedious for M > 2 
and the results are not interesting qualitatively. In this case one can simply solve the problem 
numerically. Although we derive the expressions for up to M — 3 we do not go beyond M = 2 
in the examples studied in this work. 
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Figure 2.4: Examples of the low angular momentum states tp m (0) in the short stadium. The 
effective potential V{9) is shown as a thick line. The WKB solutions of Eq. (2.10) are compared 
with the numerical solutions of Bogomolny's equation tjj = Tip (see Sec. 2.3.5). The zero axis 
for each state is its WKB effective energy E m . For the librational state the simplistic boundary 
condition ip2{0) = at the turning point was taken in the WKB case. The states and the 
potential are symmetric at 6 = 0. k^fe = 42.3 is fixed. Inset: the angular momentum 
representation of the rotational states m = 48 and m = 168 (Sec. 2.3.4). 
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Figure 2.5: The low angular momentum state n = 10, m = 2 for the short stadium with 
e = 0.05. The absolute square of the wavefunction is calculated from the perturbation theory. 
The theoretical wavefunction is inaccurate near the center of the billiard (see Sec. 2.4). The 
line segment indicates the length and position of the straight parts of the boundary. 
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Figure 2.6: Contour plot of the numerical low angular momentum state n — 25, m = 3 for the 
short stadium with e = 0.01. The parallel lines on the left indicate the length and position of 
the straight segments of the boundary. 
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In a perfect circle the change of angle after one mapping 8 — 8' = 2np/q = <d pq for a 
(p, g)-resonant orbit. (By definition, the particle makes p times around the billiard after q 
mappings.) The angular momentum l pq — (sgn p) fccos (0 p? /2) . For the perturbed circle we 
make an ansatz 

ip(6) = e u pi e + lk [ b h( )+ b2 h( )+-} . (2.18) 

The function should also have a slowly varying prefactor but we disregard it for now. We retain 
terms to order kb 2 in the phase of the T-operator and the lowest order term in the prefactor. 
The phase of T(9,9') and ip{8') can be expanded in 68' = 8' -8 + Q pq . We will see that 58' ~ b. 
Introducing the notation — d = Lq(0 P9 ) for the second derivative of the unperturbed chord 
length ( p^ ) and W pq (8) — L2{8, 8 — & pq ) we can write the integral 

J dd'T(e,d')^(e')^ J d e , ^j^A kLo ^-^ 5 ^ +kb " w ^-^ 

xe a pq (8~e pq )+ik[bf 1 (B-e pq )+bf 1 {e-e pq )se'+b 2 f2(0-e Pq )} (2.19) 

where we have taken into account that l pq = S' (Q pq ). The integral can be evaluated in the 
S$> approximation. The stationary point S8' st = bf[(6 — Q pq )/d ~ 6. 
After elementary manipulations equation ip — Tip becomes 

e ik[bf 1 (e)+b 2 f 2 (0)} = 

e ikL o (e pq )~ii pq e pq +i%+ikbf 1 (6-e pq )+ikb 2 [f 2 (8-0 pq )+w pq (e)+f[ 2 /2d\ ^ 20) 

We may try to proceed in the same fashion as in the (1,2) resonance case. However the 
assumption ip(8 — & pq ) — e^ lu ip(8) 1 in general, is not self-consistent. Indeed, it would mean 
that f[ has period Q pq in this approximation and W pq (8) + f[ 2 /2d = const. But W pq (8), in 
general, does not have period Q pq , so the latter equation cannot be satisfied. The correct 
conditions would be 

W)-h{8-Q pq ) = const, (2.21) 
h{6) - h{0 - Q pq ) = W pq {8) + f[ 2 /2d - E m (2.22) 

where E m is the constant to be determined. The first equation makes f[ p? -periodic. The 
second equation defines f% to compensate for the lack of periodicity in W pq (8). Thus it is 
essential to include f% in the wavefunction when W pq does not have period Q pq . 
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• pq 

'1 



Define a g-average of W p , 

V M = - S w ^ e + r0 ^)- ( 2 - 23 ) 

r— 1 

V q (8) has a period 0i g = 27r/q. (This follows from the 8 pr? -periodicity and the existence of 



an integer s such that sQ pq — 2irN + 6i g .) Performing the q-average on Eq. (2.22) we can 
determine 



f[(9) =±y/2d[E m -%($)]. (2.24) 
Thus f[ also has a period 2n/q. This leaves us with equation 

W) - f 2 (9 - O pq ) = W pq (0) - V q {9) = V pq (0) (2.25) 

that determines f 2 up to a Giq-periodic (or q-periodic) function, which appears in the next 
order of the perturbation theory. This equation has an explicit solution 

q-l 

Q 



1 

w) = h- h q = —J2 rV Pi( 9 - r@ pJ ( 2 - 26 ) 



where f 2q is a g-average. In the future we will, as a rule, leave out the tilde. If the perturbation 
is expanded in the Fourier series W pq (0) = J2i Wie lW then V q (9) = J2i W q ie lc > w and f 2 {9) = 
Y!i (1 - e- ae ^Y Wie m . The prime indicates that integers I divisible by q are not included 
in the sum. 

Apart from function f 2 the wavefunction ip(9) is of WKB form for the potential V q (9) : 



■ 1 1/2 

It is shown in Sec. 2.5.1 that ip(6) has a standard WKB prefactor (/{) ' in the lowest order 
(which has been added now by hand) . This form is not valid near the turning points and is not 
convenient in the classically forbidden region. Hence it would be helpful to have a differential 



equation, similar to Eq. ( 2.1 C| ) . Let 

V>(6>) = e a ^ e+ikb2f2{e) ij(e). (2.28) 

Then the new function satisfies the equation %p — Tip where the T-operator 

f(6,9') ~ \ ^ e iSo(e pq )~^se' 2 -7r] e tkb 2 [w pq (e)+f 2 (e-e pq )-f 2 (e)]^ 
' V 2iri 
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The second exponential is exp [ikb 2 V q (9)] . So, if ke | E m — V q (9) | « 1 we can repeat the 
argument in the previous section and derive a differential equation 



+ 2dk 2 e [E m - V q {9)} # = 0. 



(2.30) 



This equation is valid either when ke -C 1 or near the turning points. If ke ~ 1 and 9 is far 
from the turning points ip can be shown to satisfy this equation approximately by the direct 
substitution of the WKB form ( [2.27 ). Thus the first order wavefunction ip can always be 



determined from Eq. (2.3C). An example of an eigenfunction for the (1, 3) resonance is shown 



in Fig. 2.7. Figure 2. J depicts a two-dimensional wavefunction for the (1,4) resonance. 



2.2.4 Quantization conditions 

The eigenenergies E m are found after imposing the periodicity condition ipm{9 + 2ir) = ip m (9). 
For %l){9) it translates to 

i>(6 + 2tt) = e- i2n5 -ip{ff) (2.31) 

where <5 is the fractional part of l pq ^ 0. (For negative l pq we define [l pq ] as the closest integer 
from above and 8 = l pq — [l pq ] < 0.) If ip(9) has a WKB form we may again distinguish between 
the rotational [E m > m&xV q (9)] and librational [E m < maxV q (9)] cases. In the former case 
the condition is 

r 2v j 

(2.32) 



kb J d9^2d [E m - V q {9)\ = 2tt(to =f 8) 

where 6 is the fractional part of l pq . The double sign reflects the two possibilities in Eq. ( [2.24 ) 
In the librational case, if the tunneling between the wells can be neglected, 

1 



kb 



d9J2d [E m - V q {9)] = 7T m + 



(2.33) 



The wavefunction ip within a well is given by Eq. (2.17) up to a phase. Assuming V q (9) has q 
wells, the phase shift between the wells rp(6 + Qi q )/^p(9) — exp(— i8@x q ) [cf. Eq. ( [2.31 )]. Note 
that —S& pq /kb is the constant in Eq. ( |2.21| ). 



The second quantization condition comes from Eq. (2.20). With the help of Eqs. (2.21) 



and (2.22) we find 



Lo m {k) = kLo(Q pq ) - [lpq]Q P q + keE m + — = 2im, l pq ^ 0. 



(2.34) 
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Figure 2.7: Numerical and theoretical surface of section wavefunctions for the (1,3) resonance 
in the short stadium with k = 4032, e = 6.79 x 1CT 4 . The theoretical wavefunction is ip — 
coa(kbfi) sm(l9 + kb 2 f 2 ) ■ The fast dependence on 19 is removed by locally averaging (a) ip ainW 
and (b) tp cos 19. The effective potential Vs(9) has a nominal period 27r/3, but there are two 
symmetric wells per period for this perturbation. The angular momentum / = fccos^ is 
integer, so the function ip(9) is 27r/3-periodic (see Sec. 2.2.4). In the current example this 
function is the even combination of the second excited states in each well. 
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Figure 2.8: Contour plot of a numerically obtained state for the (1,4) resonance in the short 
stadium with k = 280.54, e = 0.01. The parallel lines on the left indicate the length and 
position of the straight segments of the boundary. 
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The energy levels E nm = fc„ m are degenerate under l pq — ► —l pq . This is a consequence of the 
time-reversal symmetry. The real degenerate eigenstates are the even and odd combinations 
of the states (|2.27 ) with ±p. Since I is large the states with ±2 do not overlap. However the 
overlap with the states from other resonances can create a transition between I and —I. We 



call it a resonance assisted tunneling. With this tunneling the states (2.27) are no longer the 
eigenstates, while their even and odd combinations are. The degeneracy between them will be 
removed. We will disregard this effect. 



2.2.5 Third order theory: ke 2 < 1 ~ ke 3/2 

If fee 3 / 2 > 1 one has to keep terms to order kb 3 in the ansatz ( [2.18| ) and T-operator. The 



additional terms in the phase of Eq. (2.19) are 
k 



ie-e pq ) 



-e Pq )6e ,2 + b 2 ti( 



b 3 f 3 



(2.35) 



Although the stationary point 



bf[{9 — Q pq )/d receives a 6 2 -order correction, it does not 



enter explicitly in the result. The above terms should be evaluated at S9' st and included in the 



phase of Eq. ( 2.2C ). The balance of the fc6 3 -terms yields 



f 3 (d) - f 3 (e - e pg ) 



^"(e M )/f , fim + dL 2 /de>] , /f/f 



(2.36) 



6d 3 d 2d 2 

where f[, f", dL^jdB' are evaluated at 6' = 6 — Q pq and c = const. Taking the g-average of 
both sides, we have an equation determining f2 qi the g-periodic part of f^, 



f-2 



^"(0 P9 )/i 2 cd , L' (Q pq ) 



2q 



6d 2 



f{ L (Q 



V a 



(2.37) 



pq ) 



where we used Eq. (2.5). This expression must have a vanishing angular average, since f^ q is 
the derivative of a periodic function, which determines c. Now /3 can be found from Eq. ( 2.36| ) 
up to a q-periodic function which is determined in the next order of the perturbation theory. 
Note that if fx is double- valued, so is f%. 



2.2.6 Non-resonant case 

The non-resonant tori that are not affected by the nearby resonances do not change their 
topology. They get perturbed proportionally to e. Unlike the resonant case, \fi does not enter 
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the classical picture. As a consequence (see Sec. 2.3.1), the odd powers of b drop out of the 



expansion (2.18). The unperturbed state has an integer angular momentum I, and after one 
mapping the angle changes by a non-resonant amount 9 — 9' = 0; = 2 (sgn I) cos -1 \l/k\. 
We start with an ansatz ip(6) = exp [U9 + ikef2(0)] and proceed as in the resonant case. 



Similar to Eq. (2.20) we find 



Jkef 2 (6) 



JkLa{ei)-il&i+i%+ike[h{6-ei)+Wi{e)} 



where Wi(9) — L 2 (9,9 — 9;). If we require 



/a(0)-/2(0-e,) = W,(0) 



we find the quantization condition 

kL (Qi) - IQi = 2tt I // 



1 



(2.38) 



(2.39) 



(2.40) 



Thus, at this order of the calculation, there is no shift in energy levels. (Note that the above 
expression is the Debye approximation to the roots of Bessel's function Ji(k) for large k and 
fixed l/k.) Expanding the perturbation in the Fourier series Wi(9) = J2 r W r e lr6 we find 

W r e* " 

T 



M«) = E r 



-ir&i 



(2.41) 



Note that W r =a = since, by definition, the perturbation has a vanishing angular average. 



Equation ( 2.41 ) would not be valid in a resonant case, when some denominators vanish, or 
if 0; is close to a rational <d pq . The same is true if W r does not drop off with r sufficiently 
fast, since exp (— ir@i) can be arbitrary close to 1 for some r. This is related to the problem 
of small denominators in classical mechanics \I% p%. If the series does not converge, the 



resonant solution is needed. Figure ^9] illustrates the wavefunction associated with the GM 
torus. 



2.3 Comments and discussion 
2.3.1 Classical interpretation 

In Sec. [O] we discussed the relation between the T-operator and the surface of section map. 
The map {V ,6') i — ► (1,9) is given by the equations V = -dS{9,9')/d9' and I = dS{9,9')/d9. 
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Figure 2.9: Numerical and theoretical surface of section wavefunctions for the non-resonant 
state on the golden mean torus. The system is a smoothed stadium with e = 1CT 4 , r\ = 2x 1CP 2 . 
For this torus 6; = ir(\/5 — 1). The factor e lW has been removed and the real part of the 
wavefunction is shown. 
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We have shown that if a curve l'(0') is mapped into 1(8) then the T-operator maps the wave- 
function ip'(8') ~ e i S l 'i e ') de ' into ip(8) = Ttp'(8) ~ e */ W<«>. The eigenstates of Bogomolny's 
equation ip = Tip are invariant under the map. Hence the curve I — li m (0) associated with 
an eigenstate is mapped on itself. We call this curve an invariant loop. Comparing with Eq. 



(2.18) we conclude that our method gives a perturbation expansion for this loop 

l inv (9) =l Q + k [bf[(9) + b 2 f' 2 (8) + ■■■] (2.42) 
and similar in the non-resonant case. Of course, this expansion can be derived by purely 



classical methods (see also Sec. 2.5.3). If we neglect terms of order kb 3 and higher, 



f inv = l ± kby/2d [E m - V q (9)] + kb 2 f 2 (6) (2.43) 

will be an approximation to the invariant loop. Under the Poincare map it will be transformed 
into a new loop l\. The area enclosed between the curves Zi nv and l\ will be proportional to 
fee 3 / 2 '. If this area is smaller than 2t:H, i.e. fee 3 / 2 ' <C 1, this approximation is good for the 
purposes of quantum mechanics. This is true even if no classical invariant loop exists. 

Clearly, li nv (9) is an intersection of an (approximate) invariant torus with the PSS. If we 
start with a point (/, 8) and propagate it under the surface of section map its images will lie on 
an invariant loop going through this point. Fig. [2.l| shows some of the invariant loops. There 
are no invariant loops in a stochastic region near the separatrix, but if its typical area is less 
than 2irh, the separatrix will be an invariant curve for a semiclassical wavepacket. The same 
is true for the secondary resonances. 

The three terms in Eq. ( 2.43| ) play different roles. Iq specifies the unperturbed torus. If this 



torus is close to a resonance a perturbation will transform its neighborhood into a resonant 
chain. The resonant islands have a size y/e and are described by the g-periodic function 
/{ oc ± \j E m — V q (0). The orbits are labeled by E m , which is therefore an approximate constant 
of the classical motion. (E m = maxy g for a separatrix.) The islands, in general, do not lie 
along a horizontal line. Instead, they form a wave of size e. This wave is given by the function 
f 2l more precisely, by the non-g-periodic part f 2 determined in the second order perturbation 
theory. The g-average f' 2q distorts the shape of the islands but does not shift them with respect 
to each other. In principle, higher order corrections can be found. In the non-resonant case 
/j = 0, and f 2 describes the distortion of the torus, which is of order e. In the special case of 
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Figure 2.10: Exact and approximate invariant loops near the (a),(b) (1, 4) and (c) (1, 28) reso- 
nances in the smoothed stadium. The exact loops (thick line) are obtained by the propagation 
of the classical map. The approximate loops (thin line) are h nv (9) that includes terms to order 
i) 3 in (a) and b 2 in (b),(c). The parameters are e = 0.05, 77 = 0.7 in (a),(b), and e — 0.01, 
T) = 0.19 in (c). 



(1,2) resonance the islands are not shifted by symmetry, so f' 2 = 0. 



We calculated the invariant loops for the smoothed stadium from Eq. (2.43). For the 



(1,2) resonance they are shown as discrete symbols in Fig. 2.2 (b). The orbits follow these 
curves closely and the deviation must be due to the higher order corrections. We can also 
formally calculate the loops for the short stadium [Fig. |2.2| (a)] although the orbits stay close 



to the loops only for a short time. Moreover, the expansion ( 2.42 ) breaks down in the higher 
orders, because V q {9) has a singular second derivative. Nevertheless, there are quantum states 



localized near these loops [Fig. 2.2 (c)]. Fig. 2.10 compares k nv (6) with the orbits near the 



(1,4) and (1,28) resonances in the smoothed stadium. Figure 2.11 shows the same for the 
non-resonant GM torus. 
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Figure 2.11: Exact and approximate golden mean torus for the smoothed stadium (parameters 
as in Fig. 2.9). The numerical classical map (dots, 20,000 iterations) is approximated by 
cos ^6/ + ef' 2 {d) (solid line). The inset enlarges a portion of the figure to demonstrate the 
degree of precision of the approximation. Presumably, this is the last torus to disappear as r\ 
is decreased. 
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2.3.2 Resonant or non-resonant? 



When q becomes large the resonant islands become small (e.g. resonance (15, 32) in Fig. 2.1). 
When the area of the islands becomes smaller then 2irh the resonant structure is disregarded 
by the quantum state. Thus, on the bigger scale, the (p, g)-chain with large q looks similar 
to a perturbed non-resonant torus. This agrees with the intuitive, as well as the number- 
theoretical |p7| notion that the irrational numbers can be approximated by rationals with 
large denominators with the increasing precision. 

Clearly, the resonant solution is needed if the librational state E m < max V q exists. This 
requires [Eq. fl2.33| )] 

? » > WTd- 

We assumed that the perturbation W pq (9), as well as V q , has vanishing angular average. As 
q increases, V q , the g-average of W pqi becomes closer to the angular average and drops off. 
Another way to say it, V q is of order of the qth Fourier coefficient W q of W pq (9) = J^i Wie . 
But Wi decreases with I, usually exponentially for analytic perturbations or as if W pq (9) has 



s — 2 continuous derivatives. If V q -C q 2 /k 2 ed then kbf[ ~ ±ky/2edE m shifts the constant part 
of the angular momentum so the nontrivial expansion begins with fce/2 as in the non-resonant 
case. 



Condition ( 2.44 ) also means that the area of a resonant island ky/ edV q /q > 1 as we implied 
above. Note that when k\fi <C 1 all phase space of the system can be treated with the non- 
resonant perturbation theory in agreement with our earlier estimates. 



The area of all resonances in the system can be estimated as S P q ^V £ dV q < 
k\JTd'Y^ q <l\fVq- So if V q drops off faster than g~ 4 (i.e. the perturbation has three continu- 
ous derivatives) the total affected area will vanish as e — > 0, as in the KAM theorem |p2|. 



The perturbation AR(6) = \sm8\ — 2/ir introduced in Sec. 2.2.1 has a discontinuous deriva- 
tive, so V q ~ q~ 2 . In this case the KAM theorem breaks down completely |24 . We expect our 
leading order solution, which does not include the derivatives, to be valid. Our theory contains 
the singular second derivatives of AR in the fce 2 -order. We assumed that fee 2 is small. How- 
ever, multiplied by an infinite V" it may bring finite corrections. We call them the diffraction 



corrections. Figure 2.4 compares the theoretical and numerically exact wavefunctions in this 
case. 
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It is easy to show that if AR(6) has only one Fourier harmonic, say AR(8) = cos (19) for 



some integer /, then V q ^ only if q is a divisor of I. (See also Sec. 2.5.2.) For example, I = 1 
and I = 2 are to first approximation a shifted circle and an integrable ellipse, respectively. 
Thus AR can be transformed to a perturbation about an integrable system with perturbation 
parameter e 2 rather than e. In the case of ellipse V2 ^ is the bouncing ball state which does 
not exist in a perfect circle. The states with the larger angular momentum are similar to the 
states in the circle. 

2.3.3 Close to a large resonance 

Another important question is whether a given unperturbed torus I (resonant or non-resonant) 
will be affected by a neighboring resonance l pq once the perturbation is on. This is the case 



if the distance \l — l pq \ < k^/edV q , the size of the resonance. We can estimate \l — l pq \ ~ 
\dl/dA8\ \&i — Q pq \ where AO is the change of angle after one mapping. Note that \dl/dA9\ = 
kd and |0j — O pq \ > q~ 2 , which is a generic precision when a number is approximated by the 
rationals [Jl7|. Thus if V q > d/eq 4 the torus I may be affected by a (p, q) resonance. If V q drops 
off faster than q~ 4 only finite number of resonances may affect the torus. 

If the torus is affected by a larger resonance, two cases are possible. If \l — l pq \ is sufficiently 
small, the topology of torus I will be changed, i.e. the torus will disappear and a new librational 
torus will appear inside the separatrix. In this case, obviously, torus I cannot be used as a 
starting point of the perturbation theory; instead, the expansion should be done near the 
resonant torus l pq . On the other hand, if \l — l pq \ is large enough, so that the topology is not 
changed, torus I still exists outside the separatrix, although it may be strongly distorted. In 
this case torus I may be considered either as a rotational trajectory for the resonance (p, q) or 
the perturbation expansion can be done directly near the torus /. In fact, the larger \l — l pq \, 
the better is the expansion near I and the worse is the expansion near l pq . Note that I can be 
another resonance (p',q') with q' 3> q. Consider, for example, the resonance (15,32) in Fig. 



2.1. It is a secondary resonance for the large (1,2) resonance. The perturbation expansion 
near the (1,2) resonance yields only fi part of the (15,32) resonance, i.e. it cannot describe 
the loops. However, for small e the loops are small and will not be reflected in the quantum 
state anyway. 
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Suppose I is a non-resonant torus. The above discussion suggests that if I lies outside of 
the separatrix of a large resonance l pq , we should be able to describe this torus by both the 
non-resonant perturbation theory for torus I and the resonant perturbation theory for torus l pq 
to some approximation. For a rotational state far enough from the resonance we can assume 



\V q \ <C E m and expand bf[ ~ ±\/2edE m (l — V q /2E m ) . We are going to show that 



l pq ± ky/2edE m (1 - V q /2E m ) + kb 2 f 2 = I + kb 2 f 2 



21 



(2.45) 



where f 2l is a non-resonant function. The constant part of /{ shifts the unperturbed angular 



momentum, so l pq ± k\J2tdE m = l pq + SI ~ I. According to Eq. (2.41) 



/»(*) = £' 



L (ej)(l + e- ir e <) AR r e* 



2(1 



-ir&i 



(2.46) 



where AR r is a Fourier component of AR(9). We can expand Eq. ( 2. 46] ) in £6 = O; — Q pq oc ^/e. 
The denominator becomes small when r is a multiple of q. We separate those terms. Then 

Lq (Op<?) 



(59 



(2.47) 



where /a(^) is the part of that vanishes under the q-average. The sum is (j)<d 1 + L /Loj V q . 
With 5Q- 1 = -kd/Sl + io"/2io + O {61) the leading term is ^V q y/d/2eE m , which takes care 



of the respective term in Eq. (2.45), when multiplied by kb 2 . The next order terms give f 2q , 
the g-average of f 2 , since, according to Eq. ( |2.37 ) , 



/2g — 



(2.48) 



_io(@pg) 2i (9 pg )_ 

We thus confirmed that a rotational torus near a large resonance can be chosen as a starting 
point for the perturbation series in the semiclassical approximation. 



2.3.4 Localization in angular momentum 



As we explained in Sec. 2.3.1 , the derivative of the phase of the eigenfunction ip{9) approximates 
the invariant loop /i nv (0) hi the PSS phase space [see Eq. ( 2.42) )] . This suggests that the angular 
momentum spectrum of a given eigenstate is concentrated within the range of / covered by its 
invariant loop. Formally, we define the wavefunction in the angular momentum representation 



as 



dO [E„ 



V a 



-V 4 if 6 d$'lhrv(.0')-l6 



(2.49) 
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If fci/e is large, the can be employed. The stationary phase condition can be satisfied 
if h nv (9) = I for some 8. In particular, for a librational state near the (1,2) resonance the 



spectrum range is < \l\ < ty'i^/max (E m — V). It has much overlap with zero angular 
momentum. In a higher resonance the range is centered near the resonant angular momentum 



l pq and spreads by ky^ey 2d max (E m — V q J in both directions. In a non-resonant state the 
spectrum is highly localized near the unperturbed angular momentum, the spread being of 
order ke. The numerical examples are shown in Fig. |2.4[ 

Outside of the range of h n v{Q)i the angular momentum components \ipi\ decay exponen- 
tially for smooth V q {9). In this case the range of l- lnv {9) defines the localization length. In the 



stadium case AR(6) ~ |sin#| (Sec. 2.2.1), on the other hand, they decay as a power law, that 
makes it necessary to define the localization length more precisely |25| . 

Note that ipi can be calculated by diagonalization of T-operator in the angular momentum 
representation. Consider the matrix elements T w = {2-k)- 1 J dddd'T(d, 9')e l ( l ' e '~ w ) . The 
stationary phase conditions require I = dS/89 and I' — —dS/d9' where the action S(9,9') is 



given by Eq. (2.3). It follows that for \l — l'\ > ke the matrix elements are small, i.e. Tw is a 



band diagonal matrix. 

2.3.5 Numerical computations 

To verify our theory we conducted several numerical checks. Essentially, we compared the 



wavefunctions and energy levels obtained in Sec. 2.2 ("theoretical" solution) with the numerical 
solution of Bogomolny's equation tp = Tip, which we presume to be semiclassically exact. The 
problem can be separated into two parts. First we solve a more general equation 

e*"» { *ty = T(Jfe)V (2.50) 

assuming k continuous (it follows from the semiclassical unitarity of T). Then we can find 
the allowed fc's from the quantization condition u m (k) = 2irn. Unless we are interested in the 



energy quantization, we may consider only the general eigenvalue problem (2.50) with fixed k. 
This is done in some numerical examples to study the approximation of the wavefunctions. 
Equation ( gjCj ) can be solved both theoretically and numerically. The theory of Sec. |2.2| 



remains unchanged but one should not equate u) m {k) to 2-7m in Eq. (2.34). To solve the 



equation numerically |3l| we choose an ansatz tp(9) and evaluate a discrete function F(n) = 
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(T n ip) (6q) at sonic 8 = 9q. We then find its Fourier image F (lu) , which consists of the peaks 
near the eigenphases u> m (k) (see below). We choose cu m of the highest peak. The eigenstate 
ifrm{0) = (( e ~ luJm T~) n ip(0)) _^ is the average over n. 

To justify this method expand the ansatz ip(9) = Yl m c m4>m{9) in the eigenstates of the T- 
operator. Then F(n) » £ m c m e^ n ip m {9 ). The Fourier transform F (ui) = F(n) e - tujn « 
(27r) _1 ^ to c m S (u> — uj m ) ip m (6o) has peaks at the eigennumbers uj m . If the ansatz is close to 
some eigenstate ijj mo (9), the respective peak will have the largest weight (for a generic 0$). 
In the combination (e-^T)" ip(6) = c mo Tp mo (6) + J2' m c rn.e z ^ m ~ UJm ^ n ip m (9) the sum over 
to ^ toq disappears after the n-average is taken. 



Figure 2.4 compares the theoretical ("WKB") and numerical wavefunctions for the short 
stadium (Sec. |2.2.1 ). The librational and rotational states near the (1,2) resonance are shown. 
The considerable difference in the bound state is due to the simplistic WKB solution of Eq. 



(2.10) with the condition that tp(9) vanishes at the turning point. This one-dimensional equa- 
tion can be solved, of course, numerically, but even the simple approach renders the main 



features of the state. Figure 2.2 (c) contains the Husimi plots of the numerical wavefunctions. 
In addition to the librational and rotational states, we show the state near the separatrix. This 
state has E m just greater than the maximum V(9), which means that the wavefunction has an 
excessive weight near an unstable periodic orbit (due to the WKB prefactor). This explains 
the observed "scars" of unstable orbits p2| . Note that scars do not appear in the previously 
developed perturbation theories that quantize the neighborhood of a stable orbit (minimum 



V). It is necessary to know the potential V near its maximum, not minimum. Figure 2.7 
shows the wavefunctions for the (1,3) resonance. In this case we remove the fast dependence 
on W by local averaging. The wavefunctions for the non-resonant GM torus in a smoothed 



stadium are shown in Fig. 2.9. The examples of the two-dimensional wave-functions are given 



in Figs. 2J5, |2J| and 2_j 



The numerical and theoretical eigenphases modulo 2ir, as defined by Eq. ( 2.50 ), are com- 



pared in Table 2.1. We consider the states with different to, but all belonging to the lowest 
(1, 2) resonance in the short stadium billiard. The wavenumber k is fixed. One can see that the 
errors in Wtheor are small compared with the average spacing, oj m +i —u m , within the resonance. 
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The spacing is of order [cf. Eqs. (p.9[ ) and ( |2.34[ )1 



ke (E m+ i - E m ) ~ y/e. 



(2.51) 



The energy levels (fc nm ) are given by the condition w m (fc) = 27m. Table 2.2 shows numer- 
ical and theoretical wavenumbers for fixed m. The error is much smaller then the differences 



7r and k 



n,m+l 



(2.52) 



given by Eq. (2.11). Note, however, that our theory gives the energy levels grouped by the 



resonances they belong to. In this example we have found the levels close to the lowest 
resonance. There are levels coming from other resonances (i.e the states with higher angular 
momentum) in the same energy range. The mean spacing of all levels in the billiard in terms 
of k is 2/k, which is of the order of the errors committed. Nevertheless, the results are still 
useful since the interaction between the states belonging to different resonances is small. 



2.3.6 Possible experiments 

A number of experimental techniques could, in principle, be used to verify the theoretical 
results (see also Sec. |4.3| ). For example, using a scanning tunnel microscope (STM) individual 
iron atoms can be positioned on a copper surface to form a boundary of a two-dimensional 
domain called a quantum corral p7| . The two-dimensional electron gas on the surface will be 
partially confined within the corral, and the spatial images of the electron density standing 
waves can be produced with STM E3 |. 

Another possible system would be a shallow container of liquid which is vibrated to produce 
standing surface waves [ij . The dissipation may limit the selection of eigenstates that can 
be observed. 

It is also feasible to utilize the analogy between a wavefunction and an electromagnetic 
wave in a resonator. Both the wavefunction in a billiard and the electric field in a cavity 
satisfy the Helmholtz equation with the Dirichlet conditions. A technique was developed |3f| 
that allows to measure the spatial distribution of the magnitude of electric field in a microwave 
cavity. 
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Table 2.1: Numerical (cj num ) and theoretical (wthoor) eigenphases modulo 2ir, compared for 
states with different m, but all belonging to the same low angular momentum resonance. 
AR corresponds to the stadium billiard. WKB "energy" parameter E m — V m - ln is also given, 
k = 1000 and e = 6.1 x 10~ 4 are fixed. 
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Table 2.2: Energies k with different quantum numbers n, but same m — 11, computed numer- 



ically and found solving Eq. (2.34). e as in Table 2.1 



2.4 Two-dimensional wavefunctions 



Although the surface of section wavefunction tp(9) together with the quantization conditions 
contain most information about the state, it might be occasionally necessary to reconstruct the 
actual two-dimensional wavefunction \1/ (r, 9) in the billiard. The procedure is straight-forward 
and involves evaluation of the integral ( 1.33| ) by the S<&. 

We start with the PSS wavefunction (or one of its branches) 

1 



u pq 8+ik[bh{8)+b 2 f 2 (e)] 



(2.53) 



where b = y/e. Without worrying about normalization, and thus leaving out constant factors, 



we write the Green function (1.34) 



G(r,6') = 



d 2 L{r,0') 



1/2 



JkL(r,6') 



(2.54) 



drj_dd' 

where L(r, 9') is the length of a direct orbit between point 9' on the boundary and point 
r = (r, 9) inside the billiard. In the unperturbed billiard there are two such orbits of angular 
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momentum l pq for any given r (Fig. 2.12). They start at the points 



a /(o) _ a ®?"z _i_ / lpq 



?i 2 = 9 - ^ ± cos- 
2 



A - ; 



(2.55) 



on the boundary. These points will satisfy the lowest order stationary phase condition in 
the integral J d9'G(r;9')ijj(9 l ). We therefore can expand the phase of the integrand in 59' = 
9' — 6'i2 ~ b. We keep terms up to b 2 . In the prefactor we keep only the lowest order. It can 
be shown that at these points 

d 2 L 7 (1) 

(2.56) 



dr ± d9' r(o) 



where 



iS=i (0) (r,< 2 ) )=7(l)T7W 
is the length of the orbits for an unperturbed circle and 

7 (r) = yfr* - (l pq /k f. 



(2.57) 



(2.58) 



Integrating over SO' near the two stationary points we find the two-dimensional wavefunction 
(up to a constant factor) 



*(r) 



E 



1^2 V^tW 



±lf , r>^, (2.59) 



ffrlsO, r<^. 



(2.60) 



Note that the theory breaks down when r — > (l pq /k) . This is easy to understand in the case 
of a perfect circle. Indeed, r = l/k is a caustic for the trajectories with angular momentum I. 
The semiclassical approximation normally fails near caustics. There are no classical trajectories 
in the region r < l/k, which means that the wavefunction should be exponentially small there. 
This is the case for the Bessel function Jiikr) for r < l/k and I 1 EJ , 

The explicit form of wavefunction Q2.59| ) is rather cumbersome in the general case. It 
simplifies significantly for the I = resonance. In this case 0'i = 9 and O 1 ^ = 9 — tt and the 
wavefunction becomes 

cosTexpiWiW (261) 
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kr + kb z 



jf _ AR (9) -AR(6- n) 
2r 2 



/ m 






4) 







(2.62) 



We assumed that u> m — irm in Eq. ( [2.11 ) and dropped all constant factors. The non-7r-periodic 
correction to the phase kb 2 f2 (9) = for the lowest resonance and the 7r-periodic kb 2 (9) can 
be included in the phase, if needed. Note that the first factor in Eq. ( |2.61 ) is the Debye 
expansion |Il| for a linear combination of Bessel's functions a (9) Ji (kr) + (3 (9) Ni (kr) with 
I = kbf[ when 

fcr>Z>l. (2.63) 



This relation will rccmcrgc in Sec. 5A where we solve the same problem by the Born-Oppcn- 
heimer approximation. Equation ( 2.63| ) sets the limits of our theory, as well. Indeed, if r is 
close to the center of the billiard, the stationary phase trajectory connecting this point with 
the boundary is not well-defined. 

It is instructive to show that the two-dimensional wavefunction ( 2.59| ) satisfies Dirichlet's 
conditions on the boundary. When we substitute r = 1 + b 2 AR (9) , we find 

, b 2 AR(9)l pq ^ b 2 AR(9)l p 



9 m 



7(7) = 7(1) + 



7(1) k : 
b 2 AR(9) 
7(1) • 



= 9-e 



7(1) 



(2.64) 
(2.65) 
(2.66) 



7(1) ' * ' W 7(1) ' 
We use these expressions in Eq. ( [2.59 ) keeping terms up to order b 2 in the phase and the lowest 
order in the prefactor. Note that the wavefunction ^ (r) tx cos [((pi — (^2) /2] , where <f>i are 
the phases of the two terms in Eq. (2.59). Using Eqs. (2.21),^] (2.22), and (2.34) one can show 
that the cosine argument is tt/2 — im, completing the proof. 

The examples of the two-dimensional wavefunctions in the short stadium are shown in Figs. 



2.5, 2. 61 and 2.8. 



2.5 Additional remarks 



"The constant in this equation is 
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2.5.1 The prefactor 



It was noted in Sec. 2 . 2 . 3| that the PSS wavefunction ip(9) has a standard WKB prefactor 
[fi (r)] ^ ^ e -^m — V^(6*)| 1. Now we lift this restriction. We solve Bogomolny's equation 



ip(9) = / d9 / T(9,9')ip(9') 



with the ansatz in general form 

i/j(9) = cxp 



U pq 9 + ikJ2b n f n (6) 



n=l 



(2.67) 



(2- 



where b = -^/e. We will be interested in the imaginary part of the phase. To evaluate the 



integral ( 2.67 ) by S$> we expand the phase of the integrand in 89' = 9' — 9 + Q pq . This phase 



is 



, 9') = kL (9, 9') + l pq 9 + kJ2 b n F n {9, 9') 



(2.69) 



n=l 



where F n (9, 9') — L n (9, 9') + /„(#). Here L n (9, 9') is the term in the expansion of the chord 
length in Eq. Q, L (9, 9') = £°° =0 b n L n (9, 9') , L n = for odd n. 



The integral ( |2.67 ) is equal to 

-1/2 



1 



\n=l 



exp i 



Qpq) k- 



(2.70) 



where all functions are evaluated at 9' = 9 + Q pq and all the derivatives are with respect to 9' . 
The prefactor can be moved to the exponent by the transformation 



i+ E & "^74' 



1 -1/2 



\n=l 



exp 



bJl 
2L» 



0(b 2 



(2.71) 



It adds the imaginary component to the phase in Eq. ( 2.70 ). This phase should be equal to 



the phase in the l.h.s. of Eq. (2.67). The imaginary part of the resulting equation is 



Im | k fj 6" [f n (9) f n (9 Q pq )} + fc 4+^3"^ } = *k + ° (&2) 



(2.72) 



We need to keep terms only to the order b. Suppose we can choose no such that kb na+1 ~ b. 
Let us make a self-consistent assumption that there is no variable imaginary part in /„ (9) for 
n < no- Collecting terms of order kb n ° +1 we can write 



Im 



fn +l(9) - /n +l( 



f f 

J 1 J ri[ 
Til 



n 



2fc6"fL / ' 



(2.73) 
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After the g-average we obtain the equation 



kb n °Imf' 



fi 



n 2f , 



(2.74) 



which means that the prefactor in function ip(9) is e~ ln V7l+°( b2 ) = l/yJJ[ + O (b 2 ) . Note 
that Im/ no is (/-periodic since, according to Eq. ( p. 72 ). Im [f no (9)— fn o {0 — ©jx?)] = 0. 



2.5.2 Fourier expanded perturbation 



Consider a system described by the perturbed Hamiltonian (2.1). Suppose the perturbation 
is expanded in Fourier series 



H 2 (I,6) = £# m (I)e i 



(2.75) 



It is known in the classical perturbation theory [[74| that each Fourier component produces a 
resonance when I satisfies the equation 



w(I) • m = 



(2.76) 



where ui(T) = 8Hq (I) /dl are the unperturbed frequencies. Below we show that the semiclas- 
sical theory agrees with this condition. 

Suppose the perturbation consists of one component H m (I) e J (™i i+™202)^ choose the 
surface of section 62 — 0, as in Sec. ^J. We show that the effective potential for the (p, q) 
resonance V q (9\) = unless LO1/0J2 = p/q = —ra^jrax. The perturbed action S2 (9,9'; E) = 
— J Hi&t, where the integral is taken over the unperturbed trajectory between two consecutive 
crossings of the PSS at 9\ — 9' and 9\ = 9. With the substitution Qy = 9' + uiyt and 92 = w^t 
we find 



S 2 (9,9';E)=- 



Tjr p im 1 6' 


p2ir 


' . ( wi 








/ exp 


1 m\ 




92 


UJ2 


Jo 


V W 2 







d9 2 



Note that H m (I) and oji (I) depend on 9 — 9' and E. The effective potential 



q-l 



V q (9) = (S 2 (9, 9 - Q pq ; E)) oc £ e ^^r P / q = 



(2.77) 



(2.78) 



r=0 



unless mip/q = I is integer. But in this case ^ exp i (jn-i^ + "^2^ $2 d&2 = 0, except 
when I = — rri2, i.e. ui ■ m = 0. In the latter case V q (9) = — ^f- H m e lmi0 . Note that the current 



definition of V q (9) differs from one in Sec. 2.2.3| by a factor of fc 



5G 



2.5.3 Direct quantization of classical perturbation theory 



In Sec. 2.3.1 we have mentioned that some results of the classical perturbation theory follow 
from our semiclassical theory. In fact, the connection between the two theories is even closer. 
In this section we are going to show how one could develop the semiclassical theory in the 
action-angle variables solely on the grounds of the classical resonant perturbation theory, 
without employing the T-operator formalism. 

First order theory 



We begin with the classical Hamiltonian (2.1) keeping terms to order e. We will quantize 



the orbits near the (p, q) resonance, for which the ratio of unperturbed frequencies loi/uj 2 = 
p/q. Following the recipe f5"2|| , we make a canonical transformation to the rotating frame 
(I, 9) — > (1, 9) defined by the equations 

h = qh, h=l2-plu (2.79) 
9 1 =q6 1 -p9 2 , 9 2 = 9 2 . (2.80) 



The new Hamiltonian is H(I, 9) = H [1(1), 9(9) j . Like the original Hamiltonian, it can be 
expanded in the perturbation series 

H(i,9)=H (I)+eH 2 (i,9). (2.81) 

The new variable 9\ is "slow" , i.e. 9\ ~ quii —pto 2 <§C 9 2 . It describes the slow deviation from 
the resonance. Under this condition one can make an infinitesimal canonical transformation 
(I, 9) — > (i, 9) — (I, 9) + O (e) that would eliminate the ^-dependence from the Hamiltonian. 
The transformed Hamiltonian has the form 

H{hA\h) = H Q (i) + eH 2 (hA\i2) (2.82) 

where the perturbed part is simply the average of H 2 over 9 2 , 

H 2 (ii,§i; I 2 ) = — — f 9 H 2 (i,9)d9 2 . (2.83) 
2tt<? Jo 

(When 8\ is fixed, the period of H 2 in variable 9 2 is 27rg.) At this stage we neglect the difference 
between the "hatted" and "barred" variables. 
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Now the problem becomes essentially one-dimensional. I 2 = ha = const is an integral of 
motion. We define ho by the resonance condition 

dH 



dh 



= 0. (2.84) 
(The meaning of this definition will become clear later, when we transform back to the original 



variables.) The energy conservation H — E requires 

^f(Ah) 2 +eH 2 {l A) = E-Eq. (2.85) 

Here Ah = h-I 10 , E Q = H (Io), and H ' = {d 2 H /dll ) 1=1q . Remarkably Eq. ( gTsg ) is similar 
to the energy conservation law for a one-dimensional particle of mass Hq with momentum 
Ah moving in a potential e_ff2(Io,#i)- We will quantize it by changing the variables to the 
appropriate operators, but first we will express the constant parameters in terms of the original 
actions Io- 

To begin with, note that 

^ = qu, 1 (l)-pw 2 (I) (2.86) 
oh 



so that condition ( 2.84 ) means that I is a resonant torus. To calculate the second derivative, 
it is convenient to express (d/dh)f 2 in terms of (d/dh) E . Using the notation from Sec. [[~7|, 
h = 9e (h) and a (I) = u>i (I) /w 2 (I) = —q'e (h) > we can write 

"\ a(I) (A.) _(') . (2 , 7) 



q\dhJi 2 \dhj l2 \dh) Il \dh 

Then the second derivative becomes 



' du 2 (qa-p) \ _ 2 , ,// 

H °- q [ dT, ) E ~ ^ W25B ll=Io (2 ' 



where we took into account that qa — p = at the resonance. If the PSS is chosen at 62 = 0, 



we can express Hq in terms of the unperturbed action (2.2). We have shown in Sec. i.7 that 
S' Q ' (Gpg) = - {2-ng%)~ 1 , therefore 

Now let us take care of the potential term in Eq. ( |2.85| ). Returning to the original variables 
we write 

B a (hJi) = ir- H 2 \I ,e 1 = —+u l t,6 2 =u 2 t \dt. (2.90) 

27rg Jo V q J 
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Note that the correction to the action between two crossings of the surface of section 
1S2 (9i,9\ — ©p<?) — — Jq^^ 2 Hidt where the integral is taken along the unperturbed orbit. 
If we define the effective potential as a g-average V q {9{) — (S2 {0\,0\ — 6 P9 )) 9 , we find that 



H 2 (i o >0i) = ~V q (6 1 = §i/q). (2.91) 

Z7T 



Equation ( |2.85 ) can be transformed to the Schrodinger equation by changing All 



id/d§i. With the help of Eqs. ( |2.89fl and ( |2.9l| ), we have 



1 d 2 ib - ~ 

+ eV 9 ((9 1 /g)V> = e£V (2.92) 



2S% (G pq ) S(0i/<Z) 2 
where 

e£ = - (2tt/w 2 ) {E - E ) . (2.93) 
It is easy to construct the complete wavefunction, which is trivial in 82, 

* {0i,92) = i> { — I e d °- 8 =$(61- -9 2 ) e a ° e . (2.94) 



q J \ q 

Setting 82 = we obtain the one-dimensional wavefunction that would also follow from the 



first order perturbation theory based on the T-operator. In fact, Eq. (2.92) is similar to Eq 



( 2.30 ) when 6\ = Oi/q. Thus both theories are equivalent to this order. 



Quantization conditions 

The quantization conditions follow from the requirement that ^ (9±, 6 2 ) be 27r-periodic in both 
variables up to the Maslov phases. With 

$ (9 + 2n) = i! (9) exp (-i2n8) , (2.95) 

where < 6 < 1 is to be determined, we obtain two equations, 

I w = l + vt+ 6, (2.96) 

Iw = n + v n --5. (2.97) 

Here / and n are integer and vi and v n are the Maslov phases. Since iio and /20 are related 



by the resonance condition (2.84), the above equations allow us to express 5, and I20 in 



terms of I and n. Note that I and n are not independent. Usually I, n 3> 5, so they are related 
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by the resonance condition on ho and ho- In the perturbed circle, for example, I is an integer 
part of l pq . 

The Schrodinger equation (2.92) with Eq. ( 2.95| ) as the boundary condition provides the 



quantized "energy" levels £ m . The total energy 



E = H (I ) - e£j^L (2.98) 

Z7T 



is found from Eq. ( [2.93 ) 



Let us compare this result with the quantization condition given by the T-operator [cf. Eq. 



S (e pq , E)-{1 + vi) Q pq + e£ m = 2ir (n + v n ) . (2.99) 

Here So (® pq , E) is the action for an unperturbed system with energy E between two crossings 
of the PSS 02 — with the difference in 9\ equal Q P q- By inverting the function So one finds 
the energy 

E = H [A0i = Q pq , S={1 + v{) Q pq + 2ir{n + v n ) - e£ m ] . (2.100) 
It can be shown that (dHo/dS) A6i — lo2/2t:. Hence expanding the above equation in e we 



obtain Eq. ( 2.98 ), taking into account that Ho (Io) = Ho [A9i = Q pq , S — Q pq ho + 2irl2o] 



Second order theory 



To find the next order correction to the wavefunction ( 2.94 ) we need to take into account the 



difference between the coordinates (I, 9) and (I, 0) defined above. In particular, A/i in Eq. 



( 2.85 ) should be changed to 

dF 

Ah = Ah - e^r- (2.101) 

001 

where the function -F(Ioi 9) can be determined from the equation fls^ 

dF - 

cj 2 - r ^H 2 (l ,9 1 )-H 2 (lo > e). (2.102) 
o9 2 

The solution of the Schrodinger equation will be modified by a factor exp (ieF) . 
The two-dimensional wavefunction becomes 

(01,0a) =i> ( j J <''"""'" ;i "" : . (2.103) 



GO 



Function F should satisfy Eq. ( 2.102| ) and have period 2tt in 8\ and 9 2 . Hence the energ 



y 



quantization does not change. Explicitly, 

F(M 2 ) = — / 2 \S 2 0i)- MkA)] d6' 2 + f 2 \e x = -) (2.104) 

^2 Jo 1 J \ q J 

where we dropped argument Io to shorten notation. The function f 2 (9±) comes from the 
T-operator theory. It is the second order term in the phase of the PSS wavefunction ip(0i) = 
ip (81) exp [UiqOi + ief 2 (81)] and it should satisfy the equation 



h (0i + & Pq ) - h = s 2 (e 1 + e pq ,9 1 ) - v q {e x ) 

H 2 {9 X ) - ffa^i,^)! . d§ 2 . (2.105) 



UJ 2 J l j a 1 =qo 1 

The 27r-periodicity in 9 X for F(8) is trivial. To show that it is 27r-periodic in 9 2 consider 

/ H 2 (§t - 2tt P , 9' 2 )de' 2 = / #2 (0i, &2 - 2ir)d6' 2 
Jo Jo 

. . . . r° 

H 2 (e 1 ,6' 2 )d9' 2 + H 2 {8 x ,Q' 2 )dQ 2 . (2.106) 

J-27T 



Then it follows from Eq. ( ^.105| ) that F(§i - 2np, 9 2 + 2tt) = F0i,6 2 ), which is the desired 
result. 

Note that F(6) is defined up to an arbitrary function -Fi(0i). This function is constant 
along the unperturbed resonant trajectory. This trajectory crosses the PSS 6 2 = at q points 
separated by angle A9 X = 2w/q. Since Fi(0i) must be single-valued on the torus (81,82) , it 
has to be 2tt -periodic in 8\. 

The above results are in agreement with the T-operator perturbation theory. To see this 
set 82 = 0. Then F (q8i,0) = f 2 (8 X ), defined up to a 27r/q , -periodic function F\ (q0\) . 

In conclusion, we have shown that the semiclassical resonant perturbation theory can be 
constructed on the basis of the classical theory without employing the T-operator. Both 
theories are equivalent, at least to the second order. The theory discussed in this section 
automatically gives the two-dimensional wavefunction. Unfortunately, the use of action-angle 
variables is rather important. The T-operator method can be more convenient in practice, 
since it is not restricted to a specific set of coordinates. Still, even here the unperturbed action 
should depend on the difference of the surface of section coordinates, which, in effect, limits 
the choice to the coordinates that are similar to action-angle. 
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2.5.4 Perturbed spherical billiard 



It is easy to generalize our theory in the case of low angular momentum resonance in a 
perturbed spherical billiard. The billiard is a three-dimensional cavity with the boundary 
r (8,<p) = 1 + eAR (9, ip) . We choose the PSS to coincide with the boundary, which is now 
two-dimensional, as well as the surface of section wavefunction ip (9, (p) . 



Following the method of Sec. 2.2.2 we expand the T-operator near 9' = n 
With the notation r = (9, ip) and r' = (9' , ip') it becomes 

8 12 V 2 . 2 



ip = ip — 7T. 



fcsin6> 
T (r, r') ~ — exp 



i2k 1 



8 



[1 + ikeV (r)] 



(2.107) 



where the effective potential V (r) = Ai? (r) + Ai? (— r) . The wavefunction is also expanded 
up to the quadratic terms. The integral J T (r, r') ip (r') d (58') d (Sip 1 ) is evaluated in the 
approximation. It should be equal to tp (r) : 



,,i'2k 



ip(-r) [1+ikeV (r)] + 



1 



1 d 2 ip 
sin 2 9 dp 2 



(2.108) 



Substitute %p (— r) 
ential equation 



e luJ ip ( r ) where u> = 2k + keEi m . Expanding in ke we obtain the differ- 

k 2 e [V (r) - Eim] ip = 0. (2.109) 



~d 2 ip 1 d 2 ip~ 
_d9 2+ sin 2 9 dip 2 _ 

Solving this equation we find possible Ei m 's. Since V (r) is symmetric the resulting wave- 
functions can be made symmetric or antisymmetric. Hence the quantization condition for k 
is 

2k + keEi m = 2%n + u>i m (2.110) 

where toi m = or ir. 

The problem of low angular momentum resonance can also be solved by the Born-Oppen- 



heimer approximation (Sec. 5.4). So far we have not succeeded in treating high angular 
momentum resonances. 



2.6 Conclusions 

In this chapter we developed the semiclassical resonant perturbation theory. The essence 
of the theory is to quantize the classical motion near the resonances. The Poincare surface 
of section method and its semiclassical analogue, Bogomolny's T-operator, effectively reduce 
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the dimensionality of the problem. The number of terms entering the perturbation series is 
determined by relationship between the parameters k and e. It would be impractical to use this 
theory if one needs to find all states in a certain energy range. Here the states are classified by 
the regions of phase space (resonances) they belong to. Each resonance is treated separately. 
We have found the expressions for the PSS wavefunctions and the energy levels within any 
given resonance. For two-dimensional systems the wavefunctions are the solutions of one- 
dimensional Schrodinger equation. So, as soon as the effective potential is known, it is easy 
to predict how the states will look like, before a numerical solution is obtained. The effective 
potential contains information about all orbits on a perturbed torus, not just the neighborhood 
of a stable orbit. This, in particular, explains scars near unstable periodic orbits. 

It should be noted that the states described by our theory are not necessarily close to the 
true eigenmodes of the system. The states associated with one resonance can mix appreciably 
with nearly degenerate states of other resonances. Although the mixing can, in principle, be 
estimated and the proper eigenmodes can be constructed, it is hard do it systematically. This 
observation, however, does not completely invalidate the theory. While the approximate states 
are not stationary, they may exist for a long time before leaking to other degenerate states. 
Such states are called the quasimodes fq]. 

There are a number of other methods that deal with perturbed integrable systems. The 
quantization of the Birkhoff-Gustavson normal form ]l5 , 38, 38| can be considered complemen- 
tary to the present method. Birkhoff-Gustavson is applied to the systems that are, essentially, 
perturbed harmonic oscillators. Our theory works for strongly anharmonic systems, as was 
remarked earlier. In addition, the Birkhoff-Gustavson method requires numerous algebraic 
transformations, which in practice makes it numerical. The adiabatic methods based on the 
Born-Oppenheimer approximation (BOA) will be discussed in Ch. [|. These methods are 
equivalent to our perturbation theory in certain cases. The BOA is rather straightforward 
and easy to understand since it is applied directly to the Schrodinger equation. It does not 
require the surface of section and automatically produces the two-dimensional wavefunction. 
Unfortunately, its application is limited to the cases where a slowly changing parameter can 
be found. For example, only the low angular momentum resonance in the perturbed circle can 
be studied by the BOA, while the perturbation theory does all of them. 
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Chapter 3 



Rectangular billiard and other systems 



In the previous chapter we developed the resonant perturbation theory for a perturbed circular 
billiard. Now we apply this theory to the perturbed rectangular billiard. In the second half 
of the chapter we consider examples of the systems without an explicit small parameter. It 
is sometimes possible to tailor the T-operator approach to these systems on the case-by-case 
basis in order to study certain classes of states, like bouncing ball or whispering gallery modes. 
The surface of section method can also be used for the scattering problems. 

3.1 Rectangular billiard — general case 

This section has mostly an auxiliary purpose. We consider a rectangular billiard with slightly 
perturbed sides. A number of particular cases of this model has been studied such as the tilted 
billiard, the long stadium, or the square in magnetic field. We will use the general results of 
this chapter when discussing these examples. The last case is particularly interesting since the 
perturbation here is the magnetic field (that breaks the time-reversal symmetry) rather than 
the shape of the boundaries. A separate chapter is devoted to this example. 

3.1.1 Classification of trajectories 

Consider a rectangle defined by its vertices (0,0), (X, 0), (X, Y), and (0, Y), where X, Y > 0. 
The boundary shape, starting with the lower boundary and going clockwise, is described by 
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Figure 3.1: Perturbed rectangular billiard. 



the functions (Fig. 3.1) 

V = e £,i(x), < x < 1, 
x = e6(y), 0<y<l, 

V = Y — e£ 3 (x) , < x < 1, 

x = X-e&fo), < v < 1. (3.1) 

In the following we assume e < 1, while ~ |^| ~ X ~ 7 ~ 1. The averages (^) = 0. Note 
that, by definition, £j > inside the unperturbed billiard. All the conclusions concerning the 
relationship between the dimensionless wavenumber k and the perturbation parameter e made 
for the circle billiard remain in force. In particular, if fee 3 / 2 <C 1, it suffices to continue the 
perturbation expansion up to the second order (assuming there are no extremely large higher 
derivatives of 

We choose the lower boundary of the billiard as the Poincare surface of section (PSS). 
Thus, an orbit might bounce off three other walls several times between the two consecutive 
mappings. We will classify the orbits according to their topology in the unperturbed billiard. 
But, first, we introduce an equivalent representation for the trajectories, using the method of 
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Figure 3.2: Method of images for the rectangular billiard, (a) tana > 2tan/3, (b) tan/3 < 
tana < 2tan/3, where a is the angle that the orbits make with the x axis, tan/3 = Y/X. The 
orbits starting within the intervals defined by points Xi and the vertices reflect from the same 
sequences of the walls. 



images. Assuming for the moment there is no perturbation, we reflect the billiard about its 
upper side and attach the image to the original. Then we reflect the resulting billiard about 
its right side to obtain a 2X x 2Y billiard. After that we continue this billiard periodically in 



x-direction. Thus, we have an infinite strip of width 2Y (Fig. 3.2). The PSS y = is identified 
with y — 2Y. All orbits go from the bottom to the top. 

Suppose an orbit makes angle a with the x-axis. We can always assume < a < ir/2, i.e. 
the ^-projection of the momentum p x > 0, in this extended scheme. The region X < x < 2X 
is responsible for the negative p x in the original billiard. Similarly, Y < y < 2Y represents 
the negative p v . With these definitions the extended unperturbed billiard has the topology of 
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a torus, and the coordinates and momenta are similar to the action-angle variables. 

We will consider only the orbits with tan a > tan (3 — Y/X, where /3 is the angle between 
the x-axis and the diagonal. In the opposite case we can switch the x- and y-directions. The 
trajectories can be divided in two classes depending on whether tana is greater (class I) or 
less (class II) than 2 tan (3. Each family of orbits with given a divides the stretch of PSS from 
to 2X into six intervals between the points < x\ < X2 < X < x^ < x^ < 2X. Orbits 
starting within the same interval reflect from the same sequence of walls as they propagate. 
Although the borders of these intervals depend on a, the sequences of walls remain the same 



within the class. Points Xi are defined in Fig. 3.2 (a) for class I and (b) for class II. Analytic 



expressions can also be written down. For example, x\ — X — 2Y/ tana for class I. The orbits 



beginning in the interval (0, x\) reflect off side 3 (as defined in Fig. 3.1) before returning to 



the PSS. For the orbits beginning in (x\ 1 X2) these are sides 3 and 4, etc. 
3.1.2 Perturbation theory 

We begin with constructing the T-operator for the system, 



T(x,x') 



A 

2TTi 



d 2 L 



dxdx' 

v(x 1 x') is the Maslov index for the orbit that goes from point x' to point x on the PSS. It 
increases by 2 on each reflection (in the original billiard). Note that each trajectory reflects off 
the upper and lower boundaries once, which changes the phase by 27r. Thus, one may count in 



v(x,x') only reflections from the side walls. As was explained in Sec. 1.6, the Maslov phase 
remains, even when we work in the extended scheme. L (x, x') is the length of the orbit. With 
the intent to calculate the two-dimensional wavefunction later, we first find C(x,y;x') , the 
length of the orbit that starts at x' on the PSS and ends at point (x, y) inside the billiard. 

We argue that the shift of the point of reflection along the boundary due to the perturbation 
can be neglected when calculating the length. Here is what we mean by that. Consider, for 



example, the right wall. It is parametrized by y [see Eq. (3.1) for wall 4]. Suppose an orbit 
with fixed end points reflects at point yo in an unperturbed billiard. If there is a perturbation, 
the orbit with the same end points will reflect at point yo + Ay where the length is extremal, 
that is dL/dy = 0. Since |^| ~ 1, the shift Ay ~ e. The correction to the length, however, is 
i (d 2 L/dy 2 ) (Ay) 2 - e 2 and can be neglected at this order of the calculation. As a result, all 
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the corrections will depend on functions & evaluated at the reflection points for unperturbed 
trajectories. 

It is convenient to express the corrections to the length C (x,y; x') in terms of auxiliary 
functions r\i (%, y; a) , i = 1, . . . , 4, where tana = y / (x — x') . To define them, consider the 
trajectory from family a that arrives at point (x, y) of the extended billiard. On its way this 
trajectory may hit the i th wall. In this case rji is equal to the value of & at the point of 
reflection. Otherwise, r\i = 0. For example, for class I trajectories, when < x < 2X and 
< y < 2Y, 



i] 2 (x,y;a) 



0, < y < a; tan a < 2Xtana, 

£2 (y — xtana) , xtana < y < Y + xtana, (3-3) 
(2Y + x tan a - y) , Y + x tan a < y < 2Y. 
iji's have period 2X in x. Note that 772 is discontinuous at wall 2, similar for other rji's. What 
complicates things even more, the extended scheme is only approximate in the perturbed 
billiard, in which the domains may slightly overlap near the boundaries, r/i is double valued 
near the i th wall where & < 0. It is equal there cither or 0, depending on which sheet the 



point is located. Equation (3.3) can still be used in the boundary regions, but the inequalities 
are incorrect there. 

Now the length becomes 

C(x,y;x') = 

\]{x - x' - 2e [772 (x, y) + 774 (x, y)}} 2 + {y - e [771 (x, y) + 2t? 3 (x, y)]} 2 = 
£o(y,a) + e£ 2 {x,y;a) (3.4) 



where 



y 

C {y,a) = - , £ 2 (a;,y;a) = -2cosa(r/ 2 + 774) - sina (771 + 2773) . (3.5) 

sin a 



When the end point is on the surface of section, i.e. substituting y = 2Y — e£i (x) = 2Y — 
erji (x, 0) , we find the length that enters the T-operator 

2Y 

L(x,x) = — h el/2 {x, a) (3-6) 

sin a 

where tana = 2Y j (x — x') and 

L 2 (x,a) = -2cosa[T] 2 (x,2Y)+r li (x,2Y)} 



G8 



- sin a [77! (x, 0) + r/x (x, 2Y) + 2 m (x, 2Y)] . (3.7) 

We will be interested in the states near a (p, q) resonance. The resonant trajectory maps 
point x' to x = x' + x pq on the upper part of the PSS, where x pq = (p/q) 2X. Hence we expand 
the phase of T-operator in Ax' = x' — x + x pq ~ i/e. By assumption, only the case when p < q 



is considered. The solution of Bogomolny's equation ip = Tip has the form [cf. Eq. (2.28)] 



ip (x) = ip (x) (*)+<?(*)] . (3.8) 



Here » sathhes the resonance condition P ,/ P , . X „/2Y or Pl . »/^ + (gf . The step- 
function g (x) compensates for the change in the Maslov phase each time the orbit crosses a 
domain boundary at x — 0, ±JC,±.2X, ... .It satisfies the condition [cf. Eq. ( 1.60| )] 



7T 

g(x) - g(x- Xpq) = --u(x,x - x pq ) (3.9) 

that has a solution 

g(x)=-ir[x/X]. nt (3.10) 
where [. . .] int is an integer part. Function ijj (x) satisfies the Schrodinger equation 

V?' + 2dk 2 e [E rn - V q (x)] ip = 0, (3.11) 

where d = (sin a) 3 /2Y, with the effective potential 

V q (x) = -{L 2 (x,a)) q . (3.12) 

The angle a has the resonant value: tana pg = qY/pX. Note that d and V q are defined with 
the opposite sign from the circular billiard. When the walls are deformed towards outside of 
the billiard, the potential tends to decrease, thus supporting a localized state. In the circular 
billiard, on the other hand, the states concentrate where the billiard is deformed inwards. 
The Schrodinger equation should be solved with the boundary condition ip {x + 2X) = 
(x) e~ l2nS , where 5 is the fractional part of 2Xp x /2n. This condition makes ip (x) periodic 
with period 2X. Once -E m 's are known, the total energy can be found from the second quan- 
tization condition 

9Vh n \2XnS\ 

- keE m = 2vrn. (3.13) 



2Yk p 
2tt- 



sma q 



2tt 
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The part of fy(x) that vanishes under (/-average is 

9-1 

<1 



l 9 " 1 r 



(3.14) 



The g-periodic part can be found from the equation similar to Eq. (2.37), but this time it will 



include the derivative of the perturbation, since Li is not of the special form (2.5) 



3.1.3 Two-dimensional wavefunction 

The wavefunction (r) = J dx'G (r, x') ip (x'), where the integral is evaluated in the stationary 
phase approximation (S$). The kernel 



G(r,x') 



1 



2V27T 



d 2 C{r-x') 



1/2 



(3.15) 



dr±dx' 

where r± is the direction perpendicular to the trajectory at point r. The Maslov phase v (r; x') 
changes by 2 on every boundary. The lowest order stationary phase condition selects one point 



x' (r) = x — y cot a p 



(3.16) 



indicating that only one classical trajectory contributes to the integral. This is the trajectory 
that approximately makes angle a pq with axis x and arrives at point r. Keeping terms to order 
e in the phase and dropping the constant factors we find the two-dimensional wavefunction 



where 



V(x,y) =ip(x,y)e t9 ( x ')- l %< r ' x ') 



<p(x,y) = ^(ir')e ip ' r+lfce {M x 0-2M £ ™-^( K 0]+ £ 2( r ; Q p9)} 



(3.17) 



(3.18) 



and x' is given by Eq. ( gig ). has period 2X in rr-direction. 

Each point in the original X x Y domain has four images in the extended billiard (up to 
a translation by 2X). Hence the physical wavefunction ^ph(x, y) is the sum of four respective 
terms ^(x,y) (Fig. |3.3| ). The combination g{x') — %v(r;x') is constant within one domain 
but differs by tt between neighboring domains. Therefore the complete wavefunction is 



% h (x, y) = <p(x, y) - <p(x, 2Y-y)- ^(2X - x,y) + <p{2X -x,2Y-y). 



(3.19) 
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x X 2X-x 2X 



Figure 3.3: Four classical orbits that contribute to the wavefunction in a rectangle. 



This function vanishes on the physical boundary of the billiard to order e. Consider, for 



example, the left border where x = e&iy)- In the first term of Eq. ( 3. 19 ) x' ~ — y cot a. In the 
third term x' ~ 2X — y cot a. The phase difference between these terms is p x [2e£ 2 (y) — 2X] + 
27T<5 — 2ke cos a [772 (0 + , y) — rj 2 (2X~,y)] , which is a multiple of 2tt. [The last term comes from 



keC 2 ] 772 (0 + ,y) = and rj 2 (2X~,y) = 0]. Thus the first and third terms in Eq. ( |3.19| ) 



cancel. The second and fourth terms cancel as well. Similarly, the wavefunction vanishes on 
the other sides of the rectangle. In particular, in order for it to vanish on the lower side, the 
quantization condition ( 3.13| ) must be satisfied. 

Setting e = 0, we find the familiar result for an unperturbed billiard: ^ p h(x,y) oc 
sin (jTrax/X) sin (irny/Y) and k mn = y (irm/X) 2 + (nn/Y) 2 . For generic X and Y the res- 
onance condition p x /p y = pX/qY is incompatible with the quantized momenta p x — ixmjX 
and p y — nn/Y for any integer p and q. This is of no surprise because the resonances are 
irrelevant for the quantum unperturbed billiard. 

3.1.4 Example: tilted square 

We illustrate the results of previous sections with the example of tilted square. This is a 
trapezoid resulted from tilting one side of the square by a small amount. This billiard was 
studied in a different context [ [fif , but no specific states were mentioned. For definitiveness, 
suppose the vertices of the trapezoid are located at points (0, 0), (0, 1), (1, 1), and (1, e). The 
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lower side y = ex, < x < 1, is chosen as the PSS. We will be interested in the states near 
(1,1) resonance, which are characterized by the diamond shaped classical orbits. The PSS 
wavefunction ifj (x) = exp (inx) "ip (x) , where k = kj \[2, indicates that the unperturbed orbits 
make angle a = 45° with the sides of the square. Note that, as in any (1, 1) resonance, function 
fi appears only in the third order theory. 



The perturbed part ip (x) satisfies the Schrodinger equation (3.11) with the effective po- 
tential in the extended scheme V(x) = \[2 \x\ (|x| < 1). V(x) is repeated with period 2. The 
requirement for tp (x) to be of period 2 makes 'tp (x + 2) = e l2K i/j (x) . Equation ( |3.11 ) is a 
piece-wise Airy equation. For the well-localized states deep inside the wells the solutions are 
approximately 

' V2, 



tp m (x) = Ai 



1/3 



< x < 1, 



(3.20) 



where £ = \/8 is the length of the unperturbed orbit, z m is an extremum (zero) of Ai(— z) for 

A *■ 1 /3 

even (odd) m and tp m (—x) = (— l) m ^ m (x) . The eigenenergy E m = z m (2£/k 2 e) and the 
total energy 



m 



1 + z„ 



2c 
nn 



2/3' 



(3.21) 



A state can be considered deep inside the well if E m <gC max V[x), i.e. z m < (V2fc 2 e/£) 



1/3 



Figure 3.4 shows the lowest surface of section state ipo (x) ■ It is localized near x = 0, 
which for the (1,1) resonance means that the two-dimensional state is concentrated along the 
diagonal y — x. To check this one looks at the two-dimensional wavefunction in the physical 
domain that has the form 

V(x,y) = e Mx+V)+ikeV2[-yE m +(y-l)V{x-y)]^ m ( x _ y j 
J re M-x-y)+ike^2[yE m + {l-y)V{x-y)\ 1 p^_ x _|_ a 
_ e i K (-x+y)+ikeV2[-yE m + (y-l)V(x+y)]^ m (_ x _ y \ 



zlK (x- y ) +l keV2[yE m + (l-y)V(x+y)]^ m ( x + y ^ Q < X,y < I. 



(3.22) 



In Figs. 3.5 and 3.6 we plotted a simplified version of this function where we neglected terms of 



order ke and set n = Trn/2. To make the resulting function satisfy the boundary conditions we 
made a substitution y — > (y — ex) /(l — ex), which does not change the order of approximation. 
The reduced function keeps the gross features of the exact function. The cross-section of the 
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Figure 3.4: Function ^o(x) for the tilted square with e = 0.03, state n = 50, m = 0. 



state along the diagonal 



^ m {2x) - ?/> m (0)cos(7mx), to even 
*(x,a;)~^ (0<x<l) (3.23) 

0, to odd 



is shown in Fig. 3.7 for a numerical wavefunction, as well as the section across another diagonal 
y = l-x 

fib, n (2x — 1), to, n both even or both odd 
(3.24) 
0, otherwise 

both numerical and theoretical. The line y — 1 — x is the diagonal of the unperturbed square, 
not the trapezoid. It crosses its boundary at x = 1/(1 + e) which explains the cut-off in the 
figure. 

Note that the neglect of the fee order terms in k leads to the symmetries ip m (— x) — 
(— l) m i/j m (x) and V'm (x + 2) = (— l) n x/) m (x) for all x S (— oo, +oo) . This is consistent with 
the fact that the potential V(x) is even at x = and x = 1 and the Schrodinger equation 



(3.11) is solved with the real boundary condition. In general, however, the complex condition 
ip (x + 2) = e l2K "ij} i x ) leaves a weaker symmetry ij) m (—x) = (— l) m ?/>^ (x) . 

We have shown that even with a small perturbation of the side of the square the strongly 
perturbed localized states can exist, as long as kyfe ^> 1. By constructing an effective potential 
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Figure 3.5: The absolute value of the theoretical wavefunction to order k^/e in the tilted square 
is plotted for the parameters of Fig. 3.4. 




Figure 3.6: Density plot for the wavefunction of Fig. 3.5. The dashed line is y = 0. 
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0.0 0.2 0.4 0.6 0.8 1.0 

x 



Figure 3.7: The cross-section of the numerical wavefunction for the state of Figs. 3.4-3.6 along 
the diagonals, (a) &(x,x), (b) ^>(x, 1 — x). The theoretical wavefunction, which is smoother, 
is also shown in (b). The diagonal of the unperturbed square y = 1 — x crosses the boundary 
at x = 1/(1 + e), which explains the cut-off in (b). 
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it is easy to see that there are also states concentrated along the side x = corresponding 
to (0, 1) resonance. We will study these states using the Born-Oppenheimer approximation 



when we revisit this billiard in Sec. 5.2. It is also clear that there are no diagonal states in the 
symmetric trapezoid with sides x = 0, 1, y = ex, 1 — ex, but there are states near x = 0. For 
the parallelogram billiard with sides x = 0, 1, y = ex, 1 + ex there are states along the long 
diagonal but no states near the edges. 

3.2 Billiards without small parameter 

In some cases the perturbation theory that we developed can be applied to a system that is 
not a small perturbation of an integrable system. However, this system may possess certain 
states that are close to the states in an integrable system. Then these states could be studied 
using the perturbation theory where the small parameter would depend on the state itself. 



As a first example, consider the ice cream cone billiard (Fig. 3.8), which is a unit circle 
for ir > \8\ > (3 and a triangular shaped region for \8\ < (3. For the surface of section states 
that are concentrated in the circular region the perturbation theory for the circle can be used. 
The perturbation of the boundary, AR(9) — 1/ cos(/3 — 9) for \0\ < (3 and zero elsewhere, 
grows as 9 goes into the triangular region. Consequently, the effective potential V q (9) is large 
near 9 = 0, repeated with period 2n/q. The PSS wavefunction obeys the Schrodinger equation 



(2.30) with e = 1. Clearly, the state will not significantly penetrate the non-circular region as 
long as E m <C max V q = q -1 sin (@ pq /2) / cos/3 [2-k jq > 2(3). This defines the small parameter 
that depends on E m . For example, the states of the low angular momentum resonance with 
small E m are concentrated near the diameter at 9 = ±7r/2. There is also an upper limit on 
the denominator q: the potential becomes shallower when 2ir/q is smaller than 2(3. 

Another example are the "bouncing ball" (BB) states in a stadium billiard || |[ |75|, [h|. 
The billiard has the straight segments of length 2a and the endcaps of radius 1, and we can 



reduce the billiard to its upper half by symmetry (Fig. 3J3). We are looking for the states that 
have low momentum parallel to the sides and large perpendicular momentum. If such states 
are concentrated outside of the endcaps, they will not be too different from the localized (0, 1)- 
states in a perturbed infinite channel or a long rectangle. The deviation of the upper boundary 
from the straight line is described by the function £(x) = 0, |a;| < a, « (|a;| — a) 2 /2, 
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Figure 3.8: Ice-cream cone billiard. 



\x\ > a. The effective potential V(x) = 2£(x) is multiplied by a large factor of k 2 in the 
Schrodinger equation, thus it is almost a square well with V — > oo for |ar| > a. Therefore 
the states are contained in the straight region if E m <C 1. Since E m ~ (m/ka) 2 , the billiard 
should be longer than m/k to apply the method. (The opposite case is a perturbed circle.) 
The ratio of the parallel and perpendicular momenta p x /Py ~ \JE m is small, as expected. 
Since we consider only half of the billiard, we need to find the states with both Dirichlet's 
and Neumann's conditions on the lower boundary. In the latter case the T-operator has an 



additional Maslov phase ±7r, which should be added to the quantization condition (3.13). The 
theoretical state n = 10, m = 1 is shown in Fig. 3.10] . The similar results for the BB states 
can be obtained by the Born-Oppenheimer approximation ||. The theory can be generalized, 
of course, to include the perturbation of the boundary of the stadium. For example, the radii 
of the endcaps can be slightly different so that the straight segments are tilted (63]. Then, 



from what we know, the BB states will be shifted towards the wider end (Fig. 3.11). When 



the sides are strongly tilted we have an ice cream cone billiard (Fig. 3.12). 

Kudrolli et al. [4S| observed the surface wave patterns in a Bunimovich stadium filled with 
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Figure 3.9: The Bunimovich stadium with the straight sides of length 2a and the endcaps of 
radius Ro = 1 reduced to the upper half by symmetry. The deviation of the upper boundary 
from the straight line is £(x). 



Figure 3.10: Square of the theoretical wavefunction for the bouncing ball state n = 10, m = 1 
in the Bunimovich stadium. 
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Figure 3.11: Bouncing ball state in the slanted stadium. The radii differ by AR/R = 0.1. The 
quantum numbers are n = 10, to = 2. 




Figure 3.12: Localized state in the "baseball" stadium. The cap radii are R\ = (1/V20K, 
R2 = (3/V20K, where I is the distance between them. If treated as a perturbed circle the 
quantum numbers are n = 10, to = 1. The theoretical wavefunction is inaccurate near the 
center of the billiard (cf. Sec. 2.4). 
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liquid. The modes were excited by oscillating the container vertically with an appropriate 
frequency. Only three types of all possible modes were actually observed, including the BB 
mode. It was argued jjj that the whispering gallery mode and the excitations of other types 
do not appear due to strong boundary dissipation. 

The notion of BB states can be generalized to the states near a short periodic orbit that 



connects opposite sides in a generic billiard |56j, |23j (Fig. 3.13). Assume constant radii of 
curvature R\ and R 2 in the neighborhoods of the points joined by the periodic orbit. Then 
the length of a nearby chord in terms of the local coordinates along the boundary is 

\2 



L\s\, s 2 ) = a h 

1 ' ' 2i?i 2R 2 



2d 



(3.25) 



where 2d is the length of the periodic orbit and s\ and s 2 measure the distance from the periodic 
orbit along the boundary with counterclockwise and clockwise positive direction, respectively 



(see Fig. 3.13 ). Let the PSS coincide with the boundary. We define the PSS wavefunctions 
locally on each side and assume that they are well localized. Then we need to solve a couple 
of Bogomolny's equations 



The formula 



ipi (si) = / T(s 1 ,s 2 )tp2 (s2)ds 2 , 
-02 (s 2 ) = J T(s 2 ,si)ijji (si)dsi. 

/OO 
e -^-y?H m (ax) dx = ^(l- a 2 ) m ' 2 H m 
-OO 



ay 



suggests looking for the solution in the form 

k 



1pi (Si) = ^/o^exp ( --^ a % s \ j Hm I \l 'j}Pi S i ) (* = !' 2 ) 



(3.26) 



(3.27) 



(3.28) 



where H m (x) is a Hermite's polynomial. The constants ai and f3i are determined after sub- 
stituting the ansatz in the Bogomolny equations. They are 



ai 



Pi = 



1 



d 


- 1 




d 


- 1 


Ri 



(±-1 

Ri J \R 2 



\ Ri 



(3.29) 
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Figure 3.13: Short periodic orbit of length 2d and a nearby orbit between points s\ and s 2 in 
a generic billiard. The positive and negative directions are indicated. 



and similar for a 2 and /3 2 - This approach is valid only if < — 1^ — lj < 1. Also the 
localization lengths y/d/kai should be smaller than the characteristic scale of the boundary, 
including iij. The quantization condition is 

kd = ir[2n+^- + - ) +(f{2m+l) (3.30) 



2 4 



where 



H^!rO a W(£- 1 )(i- 1) - 



In the case of equal radii R\ = i?2 = R we have the same functions with a — (3 = 
\Jti(2-Ti) on each side. The phase <p = | arcsin (-^ — l) . Note that when d = 2R the 
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state is dclocalized, as expected for a perfect circle. Coincidentally, we have the correct WKB 
quantization in this limit, although the theory obviously fails. In the confocal case d = 
Ri = i?2 the solution is ip (s) = exp (— ^js 2 ) H m (^J^^j with the quantization condition 
kd = tt (2n + f + |) . 



3.3 Whispering gallery modes 

A wide range of billiards can support the states that are concentrated close to the boundary 
and called the whispering gallery modes (WGMs). This effect was first discussed by G. B. Airy 
and Lord Rayleigh pi, and Keller and Rubinow BrJ obtained the quantization conditions. 



Apart from the trivial cases the WGMs should be regarded as quasimodes (see Sec. 2.6). The 
WGMs in the asymmetric resonant optical cavities with emission have a long life-time and 
might be useful in lasers . There is a recent experimental observation of the chaos-assisted 
tunneling between the WGMs in a superconducting microwave cavity p9| , as well as of the 
excitation of WGMs by a vortex in a Josephson junction jfTj. 

The standard problem can be solved by various analytical methods giving essentially the 
same leading order result. Among them the ray method by Keller and Rubinow, the parabolic 
equation and the etalon methods are reviewed in Ref. 0; the Born-Oppenheimer approxima- 



tion is discussed in Sec. 5.5; and the T-operator is used in this section. All these methods are 
based on the adiabatic assumption that the curvature of the boundary is slowly varying (the 
estimate is given below) and never vanishes. Under these conditions a classical (non-periodic) 
orbit is proven J50| to have a caustic and an adiabatic invariant, meaning that the invariant 
tori can be constructed in the perturbation theory. If a boundary has at least one point of 
zero curvature, there are trajectories infinitely close to the boundary that reverse themselves 
and thus do not have a caustic |5q| . A wider variety of classical wispering gallery trajectories 
is found in a magnetic billiard |39). It should be noted, however, that the existence of classical 
tori is a sufficient but not necessary condition for a quantum state to exist. For example, 



Fig. 1.14 shows a state localized near the boundary in a short stadium billiard. We discuss a 



possible way to deal with such systems in Sec. 3.4 



To apply the Bogomolny's equation method we need to find the length of the chord between 
two points on the boundary. The classical trajectory stays close to the boundary, so, typically, 
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Figure 3.14: Contour plot of the whispering gallery mode in the short stadium with e = 0.05 
(see Sec. 2.2.1). The numerically obtained wavefunction with k = 242.7611 is shown. There 
are no caustics in this case and almost no classical orbits that stay close to the boundary for 
a long time. The parallel lines on the left indicate the length and position of the straight 
segments of the boundary. 
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the two points will be close to each other and the boundary between them can be approximated 
by an arc of a circle. Let s be the distance along the boundary. It can be shown that the chord 
length 

u „/|3 



L(s,s') = 



s -a- - 



24 [R 

_/|3 



0(\ S ~sf) 



IS — s 



s-s'\ 4 R' (s) 



(3.32) 



2A[R(s)Y 2i[R(s)r 
where R (s) is the radius of curvature. Although the best estimate for the radius is at the 
middle-point, it is more convenient to express it locally at point s. We make an ansatz for the 
wavefunction 

V>(s) = C(s)e lk[s - f(s)] . (3.33) 

The wavefunction can be found using the first two terms in Eq. ( 3.32| ). One can show that the 
inclusion of the (small) third term in this form does not result in the R' / R-ordev correction in 
the wavefunction. Using the stationary phase approximation (S$>) in J Tipds we find 

2/3 



/'(*) = 



3^ (re- 1/4) 
AkR (s) 



and the prefactor 



C(s) oc R-^is) oc y n T{s). 



(3.34) 



(3.35) 



The quantization condition comes from the periodicity of tp (s) over the perimeter of the 
boundary C, 

k[C-f(C)+f (0)] = 2mn. (3.36) 

Quantum number m indicates the momentum along the boundary and n reflects the motion 
in the perpendicular direction. 



The last term in Eq. (3.32) is small if 



R 



< 1 



(3.37) 



where the S<& value \s — s'\ t = y/8f'R ~ n 1 / 3 R 2 / 3 k 1 / 3 is a typical hop the classical orbit 



makes near point s. Condition (3.37) becomes 

\R'\ < {kR/n) 1/3 



(3.38) 
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thus restricting the pace of variation of the curvature. 

It is convenient to describe the motion in two dimensions in the (p, s) coordinates where p is 
the distance from the boundary positive inside the billiard. These coordinates are well-defined 
close to the boundary, where p <C R(s). The two-dimensional function can be found by the 
standard technique, 

*(p,s) = / G(p,s;s')i,( S ')ds' 



(3.39) 

where G{p, s; s') oc exp [ikL (p, s; s')] . The distance L (p, s; s') between point s' on the bound- 



ary and point (p, s) inside the billiard is shown in Fig. 3.15. As in the case of a circle (Sec. 



2.4), there are two classical orbits (for fixed m and n) arriving at point (p, s) that satisfy the 



5$ condition. They leave the boundary at points 



V2R[Jf±Jf 



(3.40) 



When we sum the contributions from these two stationary points we get the wavefunction 



* 0, a) = 



(/' - * 



1/4 



COS 



2V2 



kR 



.3/2 



(3.41) 



The function is defined in the classically allowed region p < f'R. Clearly, p — f'R is an 
equation for the caustic. Note that \& (p, s) vanishes on the boundary, when p = 0. 

The p-dependent part of \& (p, s) is an asymptotic form of the Airy function |y| Ai (—z) 
for z — 2 1 / 3 (kR) 2 ^ 3 (/' — p/R) 3> 1, i.e. not too close to the caustic. Function ( 3.41 ) can be 
rewritten approximately as 



#(p,s) = (fci?r i/6 Ai 



z n exp 



i / l(s)ds/R 



(3.42) 



where I (s) = kR (1 — /') is an angular momentum about the local center of curvature and z n is 
the nth root of Ai (— z) . Alternatively, the Airy function can be expressed in terms of the Bessel 



function Ju s \ [k (R — p)] with variable index I (s) 3> 1. Hence, Eq. ( 3.42 ) can be interpreted 



locally as an eigenfunction for a circle of radius R (s) . We will return to this point in Sec. 5.5 
where the Bessel function will follow directly from the Born-Oppenheimer approximation. 



3.4 Scattering problem 

As was mentioned in the previous section, the billiards that do not support classical whispering 
gallery orbits may still have the quantum states that are localized near the boundary (perhaps 
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R(s)- p 



Figure 3.15: Classical orbit of length L{p, s; s') from point s' on the boundary to point (p, s) 
inside the billiard. The radius of curvature at point s is R(s). 



for finite energy). The candidates for this behavior are the boundaries that have short regions 
of zero, or even negative, curvature. These regions act like scatterers on the whispering gallery 
waves that can still exist along the remaining parts of the boundary. The scattering mixes 
the waves with different transverse quantum numbers n defined in the previous section. As a 
consequence, the proper stationary states (or, strictly speaking, quasimodes are the linear 
combinations of states with different n's. If only a limited number of n's contribute to an 
eigenstate, this eigenstate can be localized near the boundary, and thus look like a whispering 
gallery mode. 

We illustrate these ideas with a model of a boundary that has one point of zero curvature. 
The work is still in progress. Specifically, we assume that the curvature in the neighborhood 
of this point, say s = 0, can be described by the simplest analytic expression 

m - £ < 3 -«» 

where u is a dimensional parameter. [In the Cartesian coordinates this boundary is locally 
y(x) = (/i/4!)a: 4 where zero curvature is at x = 0.] The boundary is assumed closed in a 
regular fashion for \s\ greater than all relevant scales. 
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3.4.1 Classical theory 

Suppose a classical whispering gallery orbit is incident on the singular region from s < 0. This 
orbit generates a map in the phase space of the boundary (e, s), where e is a (small) angle 
between the orbit and the boundary at the point of reflection s. In order to derive the map 
equations, we introduce, following the recipe in Ref. Q|, the local Cartesian coordinates with 
the origin at a reflection point Si. The x axis touches the boundary such that x increases with 
s and the y axis points inside the boundary. 

Now we can write the parametric equation of the boundary 
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From this we find the angle of reflection 

y(si+i) 
x(s i+ i) 



tan a 
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[sf +1 + s 2 i+lSi + s i+1 sl ~ 3s|] + O (sV 3 



(3.45) 



where Sj + i is the next reflection point. Another equation results from the elementary geometry 

" s ' +1 ds 



(-i + i + <w 



R{s) 6 



» r s 3 



(3.46) 



Equations (3.45) and (3.46) provide the approximate classical map. The first is the cubic 
equation for s i+ i. Once s i+1 is known, it can be used in the second equation to obtain e i+1 . 
Far from the region of small curvature, where 



|Asi| = \s i+1 - Si\ < Si, 
we can expand the map equations in Asi to express 

Ae; = e i+ i - e 



V A 2 

gSiAsi 



2 As ? - 



3 Si 



(3.47) 



(3.48) 



Treating this as a differential equation we deduce that e(s) oc s 2 / 3 oc R 1 / 3 (s). This confirms 
the well-known result that 

I = eR 1/3 (3.49) 
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Figure 3.16: es~ 2 / 3 at the bounce points Sj for a classical orbit (/i = 1). The quantity is 
constant outside of the singular region near s = 0. 



is an adiabatic invariant asymptotically for e — > if R is bounded. 

In the present case, however, the adiabatic invariant is broken when the trajectory enters 



the region of small curvature (however small the initial e may be). Figure 3.16 shows the 



2/3 . . 

quantity t{S i , which is proportional to /, for a certain orbit passing through the singular 
region. Clearly, it is almost constant on the approach, then it changes significantly within a 
couple of bounces, after which it stays constant again but at a different value. This observation 



agrees with Eq. (3.47) that estimates the size of the critical region as two jumps (i.e. only one 
reflection point of any orbit lies in the critical region for this orbit). The critical region 
thus depends on the orbit. Explicitly, the size of the jump outside of the critical region is 
As(s) ~ e(s)R(s) ~ //(/is 2 ) 2 / 3 . Taking s ~ As we estimate the size of the critical region as 

|As| ~ 7 3 / 7 / M 2 / 7 . (3.50) 

This estimate will reemerge in the semiclassical regime. 

Now consider a family of orbits with the same I that are incident on the singular region 
from one side. We are interested in the "spectrum" of Fs after the orbits pass through this 



region. The final values of Ps for such family are shown in Fig. 3.17. We chose an arbitrary 
orbit with a an initial adiabatic invariant Jo, shown as a horizontal line. Then we started 100 
orbits equidistantly from the interval (sq, si), i.e. between the first two bounces of the original 
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Figure 3.17: Values of adiabatic invariant I for 100 orbits after they pass the singular region. 
The orbits had the same initial / indicated by the horizontal line. 



orbit, with the initial e = I R~ 1/3 - The fi gure shows the final I vs. the orbit number. 

Clearly, there is a considerable spread in the final values, and thus the family and the 
caustic associated with it no longer exist beyond the s = region.]] In fact, numerically, the 
maximal and minimal final values of I are related to the initial value Iq by 

Jmn/Zo « 1.96, J min // « 0.51. (3.51) 

This agrees with the relation 

Imax/Io = lo/Iwim (3.52) 

which is proven below. These results are important for the scattering theory in the next 
section. 

It should be pointed out that the right hand sides of Eq. ( [3.51 ) arc valid for any \x and any 



sufficiently small Iq. This follows from the scaling properties of the equations of motion (3.45) 



and ( 3.46 ). For any trajectory (ej, Si) there is a trajectory (e' ; , sQ = (ei, Sj/a) that satisfies the 
equations of motion in a billiard with a new parameter p! = a 3 /i. The new adiabatic invariant 

1 A generic orbit going around the generic billiard will have a certain statistical probability to increase / on 
each return to the s = region, i.e. it will make a kind of random walk in /-space. After some time it will 
eventually diffuse in the large / region, where the map equations are no longer valid. Thus there is no classical 
localization near the boundary for a generic orbit. 
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Figure 3.18: Two orbits that preserve the adiabatic invariant. 



I' = I/a 1 / 3 oc I. This proves the /^-independence. Likewise, an orbit (aej, a 1 ' 3 Sj) satisfies the 
equations of motion with the same /x, but has the adiabatic invariant a 7 / 9 /. This shows the 
independence of Iq. 

Now we can justify Eq. ( [3. 52 ). Suppose an orbit starts with Iq and ends up with I\ = 



Imnx(Io) = V^o where 77 ~ 1.96. Since there are time-reversal and s — > — s symmetries, we can 
start an orbit with I\ that ends up with Iq = rf~ I\. We claim that Iq = I m - U i(Ii)- Indeed, if 
there is a final I2 < Iq then we can start an orbit with I2 that ends up with 1\ and I\/ 12 > r\- 
But this contradicts the relation I max (I)/I — T), which is true for any / according to the 
previous paragraph. We thus showed that I m in{Ii) = V which, again, must be true for 
any I±. This completes the proof of Eq. ( 3.52) ) . 



Finally, we see from Fig. 3.17 that there are two trajectories that preserve adiabatic invari- 
ant. These are the two orbits symmetric under s — > — s. One bounces at s = 0, another at 
Si + i = —Si for some i (Fig. p. 18 ). 



3.4.2 Semiclassical theory 

The destruction of caustic prevents us from applying the standard quantization procedure, 
for example, the ray method. On the other hand, the whispering gallery mode is well-defined 
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outside the small curvature region. So, before looking for the stationary states inside a closed 
boundary, we may first consider the scattering problem of determining the amplitudes of 
transitions between the modes due to the flat region. 

Suppose a mode with a transverse quantum number n is incident from s = — oo. Its surface 
of section wavefunction ip n (s) is given by Eq. ( ^33| ). With the open boundary the energy k 



is not quantized. It is conserved throughout the motion serving as an external parameter. 



According to the Eqs. (3.37) and (3.38) the whispering gallery mode is well-defined for 



\s\ » |s-s'| st ~ (n/^k) 1 ' 1 . (3.53) 
Since the stationary phase quantity \s — s'\ st corresponds to the size of the classical jump, we 



find that the strong inequality is similar to the classical condition (3.47). In fact, the second 



estimate is associated with another classical equation (3.50). To show this, we note that, since 
the phase of ipn{s) is a projection of the semiclassical momentum onto the boundary, there is 
a relation cos e ~ 1 — f' n . Here e is the earlier defined angle between the classical orbit and the 
boundary at the reflection point s, and f' n is given by Eq. ( 3.34| ) . Then the adiabatic invariant 
can be expressed in terms of the quantum quantities as 

/^^.V^R 1 / 3 ,^ (3.54) 



where z n ss [|7r (n — j)\ 2 ^ is the nth root of Ai(— z) for large n. This connects Eqs. ( 3. 53] ) 



and (3.50). We see that the nth transverse mode is related to the classical family with the 
adiabatic invariant /„. 

As we know, the whispering gallery wavefunction ipn(s) is not valid near s = 0. It will be 
scattered into a number of transverse modes that are well-defined for sufficiently large positive 
s. We are interested in the transition amplitude between ip n (s < 0) and ip n i (s > 0). In order to 
apply the T-operator technique to the scattering problem, the surface of section can be divided 
into an external (SSE) and internal (SSI) pieces j}4|, such that ip n and ip n ' are well-defined 
on the SSE. From the previous discussion we conclude that the size of SSI is of order of two 
classical jumps. We will give an explicit definition later, but for now we assume only that the 
classical orbit that belongs to the family /„ before the scattering and to the family /„/ after 
the scattering makes one and only one bounce within the SSI. (In principle, we could extend 
the SSI to two or more bounces, which would make the theory more precise and also more 
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complicated.) 

Suppose that the SSI is an interval (— s n , s n >). The transfer operator between point s < — s r 
and point s' > s n > is 



ds"T{s',s")T(s",s) 



2ni 



d 2 L 



ds'ds 



1/2 



a ikL{s' \s) 



(3.55) 



where L(s', s) is the length of the classical orbit between s and s' that makes one bounce within 
the SSI. Note that this expression is valid only for the given n and n' in the SQ, because we 
did not include the orbits from other families that have a different number of bounces in the 
SSI. The scattering amplitude 



(T 2 ) , 

V /n't 



ds 



d s V;,( s ')T 2 ( s ', s )^(s) 



(3.56) 



When we evaluate this integral by the 5$, only those s and s' that belong to a certain 
classical orbit are selected. Namely, this is the orbit that belongs to the family I n (k) before 
the scattering and I n > (k) after the scattering. The stationary points s s % and s' st are the points 
of the last bounce before the SSI and the first bounce after the SSI, respectively. We conclude 



from Eq. (3.51) and Fig. 3.17 that, given n, the stationary points exist for those n' that satisfy 



0.51/„(fc) < /„'(&) < 1.96/„(fc), or, according to Eq. ( |3.54| ) 



Q.26z n < z n i < 3.84z„, 



(3.57) 



moreover, they come in pairs. Equation (3.57) shows that within the an incident mode 
scatters into a limited number of modes. For example, if n — 1 (the lowest mode) then 
ri = 1,...,5. 

Clearly, if the n — > n' transition is possible then the n' — > n transition is also possible 



because of the symmetry. The same conclusion follows from Eq. ( 3.57[ )f| based on the classical 
property ( 3.52| ), which is a consequence of the same symmetries. In addition, the probabilities 
of these transitions must be the same. This is ensured by the unitarity of the scattering matrix 



(3.56). In the 5*$ the T-operator is unitary and so is S„ 



The S^> result is independent of the limits of integration s n and s„< . In order to go beyond 
the 5$, one should, first of all, improve the transfer operator to have it include orbits with a 
different number of bounces. The boundaries of the SSI can be defined by the orbit in Fig. 



2 0.26 = 1/3.84 
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3.18 that does not bounce at s = 0. Let s n be the reflection point of this orbit that is the 



closest to the origin. From the map equations we find 

s n « (36/^- 2 ) 1/7 « [I087r(n - l/4)fc" V" 2 ] 1/? (3.58) 

[cf. Eq. ( 3.53] )]. With this definition the S$ orbit makes one bounce within the SSI. 

If the boundary is closed, i.e. the states ip n {s) are defined on s n i < s < —s n , then one can 
consider the propagator (A • S) n t n around the circumference. Here A„< n is a diagonal unitary 
matrix that accounts for the phase change outside of the singular region. Its elements are just 
e ifc[£„-/„(£„)+/„(o)] w h ere £ n i s the length of the regular part of the boundary for the state 
V>nO) ~ e lk ^-f^. The product (A • S) n -„ can be diagonalized. Its eigenvalues will be of 



the form e lXjtyk \ Since S n 'n ~ for n and nf outside of the range ( 3.57 ), we presume that 
the new eigenmodes will, in general, include the number of original modes within this range 
(although this qualitative statement should be verified). The wavenumber k is quantized by 
the condition Aj(fc) = 27m. So, if a billiard boundary contains a point of zero curvature, the 



whispering gallery modes may still exist (Fig. 3.14) but they are no longer simple Airy function 
type solutions. These conclusions need to be verified by the direct numerical computations, 
which arc in the plans for the future. 



3.5 Conclusions 

In this chapter we considered several examples when the perturbation theory can be applied. 
The perturbed rectangular billiard, like the perturbed circle, is one of the few systems where the 
unperturbed action depends only on the difference of their natural coordinates. It is straight- 
forward to build the perturbation theory in its general form in this case without changing to 
the action-angle variables. As we know, a perturbed system may possess the eigenstates that 
are localized near the (semiclassically) stable periodic orbits. We have an example of such 
state localized near the long diagonal in a tilted square. To have a strong localization, the 
perturbation of the boundary SL must be much greater than X 2 /L, where A is the wavelength. 

We have discussed a number of non-perturbative cases in which the T-operator method 
can be used to derive the analytic expressions for the wavefunctions and energy levels. These 
include the localized states in the ice-cream cone billiard, the bouncing ball states in a stadium 
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and near a period-2 orbit, the whispering gallery modes. Moreover, the theory allows to include 
an additional perturbation to such systems, for example, to tilt the sides of the stadium. 

The standard whispering gallery modes result from the quantization of classical families of 
orbits that are localized near the boundary. The families lie on the invariant tori that exist for 
the smooth boundaries with positive curvature in the limit of small sliding angle. When the 
boundary contains a point or a short interval of zero curvature, the classical tori do not exist 
but the localized quantum states are still possible. The bound state problem involves finding 
the scattering amplitudes, which is another area where the semiclassical surface of section 
method can be used. 
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Chapter 4 



Square billiard in magnetic field 



In this chapter we consider a square billiard where the perturbation is caused not by the bound- 
ary, but by the magnetic flux perpendicular to its plane. We concentrate on two particular 
cases: a uniform magnetic field and an Aharonov-Bohm flux line (ABFL). The latter system 
is purely quantal — the famous Aharonov-Bohm effect Q results from the phase interference 
and has no classical analogue. We will be interested precisely in the effect of phase interference 
on the cigcnstates in the billiard, and therefore even in the uniform field case will neglect the 
Lorentz force (below we estimate the field where it can be done) . We know from the previous 
chapters that a relatively small non-integrable perturbation of an integrable square makes it 



possible for the strongly localized states to exist (Figs. 3.5, p.6[ ). This is still the case when 
the perturbation is the magnetic flux J57J, and these states will be the primary focus of our 
discussion. 

Although we confine our attention to the weakly non-integrable billiard, the ergodic be- 
havior of the system for larger fields can be of interest. As in the Sinai billiard, the magnetic 
field imposes its circular symmetry on the square, which leads to chaos. Classically, depending 
on the field, mixed and chaotic regimes are possible Quantally, one is interested in the 
energy level statistics. The energy level distribution in a strongly chaotic billiard depends on 
the symmetries of the system, in particular, the time-reversal symmetry. The magnetic field 
breaks this symmetry, and, in general, the statistics will be different from the non-magnetic 
chaotic billiard Jl2[ . (The square in uniform field still has an antiunitary symmetry, however.) 
The spectral statistics of a rectangular billiard with the ABFL at the center is discussed in 
Refs. p8l|6l. 
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The measurements of the orbital magnetic susceptibility in the array of GaAs ballistic 
squares |5lJ motivated a certain amount of theoretical activity in the last few years |67| . The 
trace formula is a powerful semiclassical approach that allows one to express the density of 
states (and thus the susceptibility) in terms of the actions of the periodic orbits. For the 
relevant temperatures only the short orbits need to be taken into account. In the weak field 
regime, when the perturbed Berry- Tabor trace formula can be used, our perturbation theory 
also allows to find the susceptibility. 

The main reason to use the perturbation theory, however, is to predict and give the ap- 
proximate expressions for the wavefunctions. Even the relatively simple analytical formulas 
produce the states with the non-trivial probability and current distributions. It is possible to 
have current loops of opposite sense within one state. For a given field one finds the states 
with the overall current ranging from paramagnetic to diamagnetic, including those with the 
small overall current but strong local currents nearly canceling each other. 

4.1 Uniform field 

4.1.1 Limits of applicability 

A free particle of charge e moves in the uniform field B along a circular orbit of cyclotron 
radius R c — mcv/eB where v is its speed. If confined by the square boundary the orbit will 
consist of the respective arcs. We are interested in the regime when the Lorentz force can be 
neglected in the lowest order and the orbits are made up of nearly straight segments. This 
means that e = L/ i? c , where L is the side of the square, is a small parameter. In this case the 
classical momentum mv » eA/c, where A is the vector potential, so one can approximate 
mv ~ hk in the semiclassical picture. Then we can express e = 2tt 4> / 'tfiakL, where = BL 2 
is the flux through the billiard and (f>o = hc/e is the flux quantum. We shall see that e is the 
small parameter of the perturbation theory. 

In the units fi = c= e = L = lwe have e = B/k <§; 1 where the dimensionless field B = <f) = 
2-K<f>/4>Q- Without loss of generality we assume B > 0. In the semiclassical approximation we 
take the wavenumber k 3> 1. As explained in the beginning of Ch. |^, the resonant perturbation 
theory is required when k^fe = \JkB > 1, i.e. when the magnetic field is relatively large. The 
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upper boundary for B is given by the number of terms kept in the perturbation series and by 
the smoothness of the perturbation (which determines the onset of diffraction). If the terms 
of order no higher than M are kept then kb M+1 -C 1 where b = s/e. To keep the diffraction 
effects small we limit kb A <C 1 (see the next section), but if we do not go higher than the 
second order then the stricter condition kb 3 <C 1 must be imposed. Note that if k is fixed then 
B may not be too big, if B is fixed then k should be large enough, and if e is fixed then the 
energy is bounded. With the ABFL the diffraction is stronger, and kb 2 <C I in that case. 



4.1.2 Effective potential and eigenstates 

As far as the formal theory is concerned, the square in magnetic field is similar to the square 
with the perturbed boundary (Sec. 3T).[] Again we use the method of images reflecting the 



square about its boundaries. The magnetic field has opposite signs in the adjacent squares 



and the orbits are curved differently (Fig. 4.1). Note that we chose the origin at the center 



of the square. The Poincare surface of section (PSS) is the bottom side y = —1/2 identified 
with y = 3/2. We will be interested in the states near resonances (periodic orbits) since they 
differ the most from the states of the unperturbed billiard. Figure |4.1| , for example, shows a 
trajectory near the (1,1) resonance. This trajectory is paramagnetic, the square boundaries 
make it circulate in the opposite direction from a free particle in the field. 
The action along the orbit in the extended scheme is 



S(x,x') = k\jA+ {x -x'f + ®(x,x') (4.1) 

where $ = (e/c) J A • dr. The curvature of the orbit can be neglected because it makes a 
correction of e 2 . We will primarily consider the (1, 1) resonance since, as will become clear, it 
is the most important. In this case the PSS wavefunction has the form ip (x) = e %Kx ip m (x) 
where k = k cos 45° = k/V2. The slower varying function ip m (x) satisfies the Schrodinger 
equation 

^ + [E m - V (x)] 4> m = 0. (4.2) 

The effective potential V (x) — — fc$ (x, x — 2) / C is proportional to the perturbed part of the 
action evaluated along the unperturbed orbit of length C = V8. Note that $ (x, x — 2) is just 

lr The same lowest order results are given by the channeling method (Sec. B.3I). 
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Figure 4.1: Method of images for the square in magnetic field. The magnetic field has different 
signs in the adjacent squares and the orbits are curved differently. R c is the cyclotron radius. 
The curvature is exaggerated for clarity. The Poincare surface of section is y = — 1/2 identified 
with y = 3/2. The orbit shown is close to the (1,1) periodic orbit. It goes from x' to x + 2 
identified with point x in the original square. The square causes the particle to circulate in 
the paramagnetic sense, which is opposite to the direction of a free orbit in the field. 
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the flux enclosed by the unperturbed periodic orbit in the original square, equal B times the 
area of the loop with the positive sign for counter-clockwise orbits. The potential is then 

f -Bk (i - 2x 2 ) /£, xe [-ii] 
V{x)= { V2 ' 1 2 21 (4.3) 

! Bk [±-2(x + l) 2 j /£, [|, §] 

with the periodicity V (x + 2) = V" (x) . The potential is proportional to k 2 e, which justifies 
the use of e as a small parameter. V (x) consists of alternating parabolic wells and barriers of 
height Bk/2 joined smoothly at x = ±1/2, where the potential is odd. Note that the potential 



for the tilted square of Sec. 3.1.4 for the same resonance is even on the boundary because of 
the time-reversal symmetry. 

The potential V (x) has a discontinuous second derivative at x = ±1/2 due to the diffraction 
from the corners. It will bring a step-function to the fc6 4 -order of the perturbation theory. Since 
the perturbation expansion takes place in the phase of Bogomolny's equation, the discontinuity 
allows to shift the pieces of V" (x) independently by integer number of 2ir. This is the reason 
for the condition kb A -C 1. (By the same token, we must require kb 3 <C 1 in the tilted square.) 

The potential is negative for the trajectories leaving the PSS with the positive velocity 
projection v x . These trajectories go counter-clockwise around the square. The minimum at 
x = indicates a stable periodic orbit. In its neighborhood the potential is of the harmonic 
oscillator type, so for the low-lying states the energies are approximately given by 



lBk ( IN / 8Bk 
Em =-2T + { m+ 2 V— ' (44) 



The formula holds for m <C \J Bk/C. The lowest wavefunction tfjQ (x) = e V 2c x is localized 



near x = 0. The two-dimensional state is localized near the stable periodic orbit (Fig. 4.4). Not 
surprisingly, it maintains the paramagnetic current, as explained below. The wavefunctions 
with higher ra's penetrate into the region 1/2 < |x| < 1 more and more. Hence these states 
include the diamagnetic orbits. The states with large m are diamagnetic. The maximum at 
x = 1 marks an unstable periodic orbit, which is almost a time- reversed stable orbit. Our 



method accounts for the scars of unstable orbits (Sec. p.3.5| ), i.e. the states with E m w V r 



max 



(Fig. 4.11). 



In general, Eq. (4.2) should be solved with the periodic condition ip (x + 2) = ijj (x) . It is 
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accompanied by the quantization 

k nm = 2im/£ + E m (k) jk (4.5) 

that follows from Bogomolny's equation. Thus n is roughly the number of wavelengths along 
the periodic orbit. In the first approximation k in the r.h.s. can be replaced by 2im/£. Note 



that E m /k weakly depends on n [cf. Eq. (4.4)]. The energy 



k 2 = (2tto/£) 2 + 2E m (4.6) 



depends on B only via E v 



The two-dimensional wavefunction could be found similarly to the tilted square (Sec. 3.1.4) 



But here we do it by another method that takes into account the symmetries of the problem. 



The PSS wavefunction ip (x) = e %Kx ip m (x) can be propagated with the Green's function ( 3.15 ) 
to give a two-dimensional wavefunction 

*o (x, y) = e^ x+ y^ m (x - y - (4.7) 

where we ignored the terms of order ke in the phase. In this approximation k = im/2 
and ipm (x) solves the Schrodinger equation with the boundary condition tp m (x + 2) = 
( — 1)" tp m (x) . By symmetry, there are three more states with the same energy obtained from 
by 90° rotations. For example, ^i(x,y) = lZ^${x,y) = 1 $ Q (y,—x) where TZ : (x,y) — * 
(y, —x) is a rotation operator. The gauge is chosen to make the Hamiltonian invariant under 
TZ. The eigenstates in the physical domain will be the linear combinations of 'Si's that are also 
the eigenvectors of TZ with the eigenvalues i~ r , r = 0, . . . , 3. The eigenfunction of symmetry r 
is 



•''<••,„(,• // -,)■ (4-8) 

The eigenstate symmetry r depends on the quantum numbers n and m. By the boundary 
condition \&( r j (x, —1/2) has to vanish. The s = 0,3 terms are proportional to e iKX and must 
cancel each other. This implies ip m (—x) = —i~ 3r e~ lK ip m (x) . Since the potential V (x) is even 
and ip m satisfies real boundary conditions we can choose ip m (—x) = (—1)™ ip m (x) . Then n is 
related to m and r as 

nmod4 = [2 (1 - to mod 2) + r] mod 4. (4.9) 
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For fixed m the successive values of n cycle through the representations of 1Z. Note that 
-0 m (x + 2) = (-l) r 4> m (x) . 

These results are easy to generalize to the other resonances (p, q) . In the weak field ap- 
proximation the orbits with even pq do not enclose any flux. These orbits can be subdivided 
into pq smaller loops, half of which are positive, and half are negative, so the total flux is 
zero. (If the curvature is taken into account, there will be a small flux). In particular, the 
bouncing ball (0, 1) orbits are non-magnetic, and (1, 1) is the lowest resonance that responds 
to the magnetic field. If pq is odd, the effective potential is 

q 



Vpq (X) 

P 
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q(x-\ 




1" 


V 




4) 










to 



(4.10) 



where C pq — 2y / p 2 + q 2 and V is given by Eq. (4.3). The potential has a period 2/q. If kB <C 



pq£p q , the potential does not support localized states and can be ignored. The PSS states 
should satisfy the boundary condition %j) m {x + 2) = e~ l27r<5 7/> m (x) where 6 is the fractional part 
of np/ (p 2 + q 2 ) . The energy is 

k 2 nm « {2*n/C m f + (l + ^) E$ (4.11) 

where we have to distinguish between E$ for the resonance (p,q) and Em = (p/l) 2 Em 
for the resonance (q,p) ■ These resonances produce nearly uncoupled states with the same 
energy that are related by a 90° rotation. The splitting of the levels is negligible unless these 
resonances are close in the phase space. The latter occurs when p/q is close to 1 so that the 
tunneling can be assisted by the (1,1) resonance, but not too close to have (p,q) and (q,p) 
resonances destroyed by the (1, 1) resonance. 

4.1.3 Magnetic response 

The magnetic susceptibility at high temperatures is dominated by the Landau term coming 
from the smooth part of the density of states. At low temperatures the oscillatory corrections 
become large as the Fermi wavenumber kp increases. This contribution to the susceptibility 
can be found using the trace formula [p7| , which does not require the knowledge of the wave- 
functions or the spectrum. When the magnetic field is weak enough, the resonant perturbation 
theory can be used to calculate the share of a particular resonance in the susceptibility. This 



method is related to the perturbed Berry- Tabor formula (Sec. 3.4) 
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We start with a grand potential for a system of non-interacting electrons 

ft (T, /i, B) = -k B T In [l + e-^»-' 1 )/ feBT (4.12) 

a 

where fj, is the chemical potential and the sum is over the single-electron states. If we are 
interested in the part of the susceptibility coming from the (1,1) resonance, we include in the 



sum the subset of states given by Eq. (4.6). The magnetization M. = —dVt/dB is 

M(T, f ,,B) = -Y / ^^fD[Enm(B))- (4-13) 

nm 

Here fo is the Fermi-Dirac distribution function. The sum over n can be done using the 
Poisson sum formula. The result is 



a. 



_ . 1 sirih (its CksT/ 2k f) 



Cs\ k F - 

kF 



(4.14) 



where a m — 8E m /dB evaluated at k — kp. For the other resonances (p,q) one should sub- 
stitute E m — ► ^1 + Em 1 /2 if is odd; the resonances with even are non-magnetic. 
The magnetization decreases exponentially as C becomes large. Therefore the (1,1) resonance 
gives the strongest magnetic response. If we neglect the higher resonances, we do not need 
to keep the terms with s > 1 either. In the trace formula approach s gives the number of 
repetitions of the primitive periodic orbit. The susceptibility is given by x = dM/dB. It is 
exponentially suppressed for high temperature giving place to the Landau susceptibility that 
decays as a power law. Since x is proportional k^J 2 {kp coming from a rn and of order \fkp 
terms in the sum over m) it is greater than the kp independent Landau term at low tempera- 
ture, which is of order of one in our units. In the experiment |5l| the temperature was about 
10 level spacings, i.e. ksT ~ 20ir. The squares contained 2 6 x 10 4 electrons which makes 
kp ~ 300 -r 600. Then TrsCksT '/2kp ~ 0.5 -j- 1, and the exponential suppression is not too 
strong. If an experiment is conducted on an array of squares the result should be averaged 
over the distribution of sizes. Because of the oscillating sine factor the result will be reduced. 
Reference |67| provides the corrections. 

The magnetization depends on the response of individual states a m . Note that dE m /dB — 
d(H) m /dB = U m V $ m ) where B = Bk/2£, V (x) = V (x) /B, and H is the effective 
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Hamiltonian in the Schrodinger equation (4.2). Semiclassically, 

-1/2 



dE m 



J dxV (x) E-V (x) 



-1/2 



(4.15) 



Jdx E-V(x) 

If E = E m /B is treated as a continuous variable, the values of dE m /dB as a function of E 



will fall on a continuous curve independent of B (Fig. 4.2). The integration is between the 
turning points. The states with low m's have negative a m and thus are paramagnetic. The 
wavefunction has more weight in the region x £ ["515] which is classically identified with the 
counter-clockwise orbits and where V is negative. The states with higher m's are diamagnetic. 
In this case the weight shifts towards the turning points located in the region \x\ S [5, l] where 
V is positive. If E = V ma ,x, the curve has a cusp. The simple WKB formula fails in this case. 
The wavefunction is large near the unstable periodic orbit, which is strongly diamagnetic. As 
m gets even bigger, the wavefunction becomes more uniformly distributed in x with slightly 
higher weight near x = 1 . So the states remain diamagnetic but their magnetization goes down. 



The quantized transverse "energy" is approximately given by the Bohr-Sommerfeld condi- 



tion 



dxy/E m - V (z) = 7T ( m + - ) I V B. 



1 



(4.16) 

The number of transverse modes in the resonant island, that can be estimated as 0.6\/b~, 



depends on B. In the example of Fig. 4.2 k 142 and B — 25, and there are 15 states inside 
the well. With these parameters we have two special states: m = 14 which is close to the 
diamagnetic maximum and m — 10 which is almost non-magnetic. 



4.1.4 Visual representation of eigenstates 

As we learned in the previous sections, the square in magnetic field can support the eigenstates 
with different degree of localization and a wide range of magnetic properties. Here we show 
some examples of the wavefunctions and current distributions calculated using the perturbation 
theory (referred to as "theoretical" results) and compared with the numerical solution of the 
Schrodinger equation in the square ("numerical" results). We restricted our attention to the 
states of the (1, 1) resonance that are invariant under the 90° rotation, i.e. having the symmetry 
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1 .0 




S 



Figure 4.2: dE m /dB as a function of E. The stars are at the (1, 1) resonance states n = 62, 
to = 0, 2, ... , m max = 14, . . . Magnetic field B = 25. The continuous curve is given by Eq. 
(4.15). 



r = [Eq. fl4.8Q ]. These are the states with either even n/2 and odd m or odd n/2 and even 
m . 

A couple of comments on the calculation procedures are in order. In the perturbation theory 



the one-dimensional Schrodinger equation (4.2) was solved numerically. The two-dimensional 
wavefunction was constructed using Eq. ( |4.8[ ). We call these results "theoretical." 

In the "numerical" calculations the two-dimensional Schrodinger equation was solved by 
diagonalization of the Hamiltonian in the basis of states in the square without magnetic field. 
The states with the required symmetry are 



y/2 [cos npx cos nqy + cos irqx cos Trpy] , p ^ q odd 

2 cos irpx cos irqy, p — q odd 

\pl [sin npx sin nqy — sin nqx sin npy] , p ^ q even 



(4.17) 



The labels p, q should not be confused with the resonances. Expand ip m (x) = $ m ,l e x i n 
Eq. (4.S) where I is integer for r even and half-odd integer for r odd. Also tp m .i = (— l) m ip m ,-i- 
We expect the largest weights in the state (n, m) come from the B — states with p = § + I, 
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Figure 4.3: Decomposition of the state m = 14, n = 62, B = 25 into B = basis states. 
The area of a circle is proportional to the square of the amplitude. The numerically exact 
results are shown below the diagonal p = q. Above the diagonal are the weights of theoretical 
wavefunction Eq. (4.8). A part of the constant energy circle p 2 + q 2 = const is shown as well. 



q = 2 — I. Figure 4.3 shows these states aligned along the diagonal in (p, q) plane. The energy 
of the basis states is p 2 + q 2 = + 2l 2 dropping a factor of 7r 2 . If I 2 <C n, it is much closer 
to the base energy n 2 /2 than the next base energy [n ± 2) 2 /2 w %- ± 2n. Nevertheless, there 



are basis states with other base energies lying near the constant energy circle (Fig. 4.3), for 
example, the states with p' = + I', q' = ,J 4^ — V where I' 2 ~ n. However, the matrix 
elements of the Hamiltonian TL pq .p> q > are small if the differences \p — p'\, \p — q'\, etc., are 
large. When the condition n^> I 2 ~ m 2 breaks down, the basis states with other base energies 
will be contributing to the eigenstate. This is the case when the terms left out in the phase 
of Eq. ( |4.7| ) become important. The terms that shift n are of order {E m — Vmin) /k ~ m 2 /n. 
Thus, for large m our perturbation theory is not too accurate. 
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The current density is given by 

j (as, y) = 2Re** (x, y) [-iV - A (x, y)] * (x, y) . (4.18) 

The property V • j = ensures that the streamlines are the closed loops. However, because 
of the approximate nature of both numerical and theoretical calculations, the small left-over 
divergence in the regions of small current density may have a big effect. This makes it problem- 
atic to construct a reasonable approximation for the streamlines. To overcome this difficulty 
we imposed the divergence-free character of the current by representing j (x, y) = V X z% (x, y) 
where the function \ {x, y) = J^ 2 j x (x, y') dy' . Then the streamlines of j are the contour lines 
of x- In the exact problem x ( x i y) — X (Vi x ) ■ To see this, note that the Hamiltonian has a 
symmetry Ti (x, y) = TL* (y, x) , which means that for a non-degenerate state the wavefunction 
has the same symmetry (up to a phase factor) . From this follows the symmetry for the current 
jx (x, y) = -j y (y, x) . Then x (y, x) = - f* /2 j y (x' , y) dx' = f* /2 dx' = x (x, y) , using 

the condition j x \\,y') = 0. In practice, however, x 1S n °t symmetric. This is the pay-off for 
making the approximate current non-divergent. Therefore in our calculations we used the sym- 
metrized version | [x (x, y) + X (y> x )] ■ In the figures below we have a good agreement between 
thus obtained streamlines and other types of current density representations. 

Qualitatively, the (1, 1) resonance states can be divided in four classes. The first class con- 
sists of the states of the low transverse modes with E m — V m i n <C Bk. These states are strongly 
localized along the stable orbit and have relatively simple paramagnetic current patterns. In 



Figs. [4.4| - [4.8| the magnetic field is B — 31.4. Figure |4.4| shows a three-dimensional representa- 
tion of the numerical wavefunction n = 62, m — 0. The theoretically estimated wavenumber 
is k pa 27r62/£ w 140. The parameter k^fe = V Bk w 66 is large indicating the strong localiza- 



tion due to the resonance. Because of the localization only one of the four terms in Eq. (4.8) 



dominates each side of the square periodic orbit. For example, near x = —y = 1/4 only the 



s = term is appreciable. Here | SI/ (x, y)\ 2 ~ ipo (x — y — i) 



which is well approximated by 



a Gaussian. Near the square boundary there is interference between two terms. For example, 
near y = —1/2 the s = 0, 3 terms are large, so 



l*(z,y)| 2 



* ( ™ /2) ^o (x - y - + e-*<™/2)^ f- x - y - 
i> (x)cos[—y) . (4.19) 
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Figure 4.4: Absolute square of the numerical wavefunction for the localized paramagnetic state 
m = 0, n = 62, B = 31.4. According to Eq. (4.19), the profile along the side of the billiard is 
proportional to \ipo\ 2 and there are fast interference oscillations in the transverse direction. 



The cosine factor produces interference fringes parallel to the border. At the first maximum 
1$ (0, — g + i) | 2 is four times larger than |\& (|, — ~) | 2 . The cross-section |\& (a;, — | + -) | 2 
is also a Gaussian. The interference pattern can be seen in the density plot in the lower left 



quarter of Fig. 4.5. The rest of the figure shows the current for this state. Starting on the lower 
left and going counter-clockwise we show a density plot of l^l 2 , the streamlines, the absolute 
value of the current density |j (x, y)\ , and the vector field j (x, y) (the dot is at the calculated 
point and the stick size and direction represent the current density magnitude and direction). 
Of course, close to the boundary the current is parallel to the boundary. Its interference fringes 
mimic those of the wavefunction. 



Figures 4.6-4.8 show another state of the same class n = 60, m = 1. Its surface of section 
wavefunction is close to the first excited state of a harmonic oscillator. Clearly, the probability 
and current densities are very small near the point x = —y = 1/4. They are not exactly 
zero, however, because the terms with s = 1-3 in Eq. (14.81) have exponentially small tails 



there. Figures 4.6, 4.7 were obtained numerically and Fig. 4.S shows the results of theoretical 
calculations. We see that the theory works well in the low m case. 
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Figure 4.5: Numerical wavefunction and current for the state in Fig. 4.4. The state is sym- 
metric under the 90° rotation. Counter-clockwise from the lower left: a density plot of l^l 2 , 
current streamlines, a density plot of |j|, and a vector field representation of the current (the 
dot is at the calculated point, the size and direction of the stick represents the current density 
magnitude and direction). 
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Figure 4.6: Absolute square of the numerical wavefunction for the state m = 1, n = 60, 
B = 31.4. 
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Figure 4.7: Numerical wavefunction and current for the state in Fig. 4.6. The same represen- 
tation as Fig. 4.5. 
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Figure 4.8: Theoretical wavefunction and current for the state in Fig. 4.6. The same repre- 
sentation as Fig. 4.5. 



Ill 



The states with \E m \ <C V max = Bk/2£ belong to the second class. These states are de- 
localized and have rather complicated probability and current distributions because all four 



terms in Eq. (4.8) are important throughout the square. These states lie on the border of para- 
and diamagnetism and are weakly magnetic. Figure 4.9 shows such state n = 62, m = 10, 



B = 25. It has very small magnetization (Fig. 4.2) because the loops of rather strong dia- 
magnetic and paramagnetic currents nearly cancel each other. The figure shows the numerical 
results as before, except for the upper left corner, where the theoretical streamlines are shown. 
The theory and numerics differ primarily in the lower current regions. To understand the 
complexity of the higher m states and the lack of localization consider the normal derivative 



62,w /dy\y=_ 1 / 2 on t ne boundary (Fig. 4.10). In the low m state it would be proportional 
to \ip m (x) | [cf. Eq. (4.19)1- When m becomes larger ip m (x) extends significantly beyond the 



physical domain [— |, |1 (Fig. 4.1C ), so it has to be "folded" back to this domain with an 
additional rapidly varying phase factor e lKX . Formally, it means that the terms with s = 1,2 



in Eq. (4.8) cannot be neglected. 

The third class includes the diamagnetic states with E m near the top of the potential. The 
state at the maximum of the magnetization curve (Fig. |4.2| ) with n — 62, m = 14, B = 25 
is shown in Figs. 4.11 (numerics) and 4. 1 2| (theory). In the case of large m the theory does 
relatively poorly, although it is still qualitatively correct. The unperturbed basis states beyond 



the diagonal p + q = n appear in the Hamiltonian matrix (Fig. 4.3, below the line p — q), while 



the theory does not reflect this (Fig. 4.2, above the line p = q). This state is localized near the 
unstable orbit, i.e. this is a scar state that we predict. One can even see the two channels in 



Fig. 4.11 because if E m is just below V max the wavefunction ip m (x) has the excessive weight 
near the turning points on both sides of the maximum of the potential at x = 1. In this case 
we find E m = 628.57 and Vm ax = 648.64. For comparison, we show the state n = 62, m = 14, 
B = 31.4 (Fig. 4.13 ). Here E m = 703.4 is much lower than Vmax = 792.2, so the two channels 
are further apart. 

One can observe the change of the topology of the current from paramagnetic to diamag- 



netic by looking at the sequence of states m = 6, 8, 10, 12, n — 62, and B = 25 (Fig. 4.14 ). 
Similar effect could be achieved by the change of magnetic field. 

The states with E m ^> Vi nax belong to the fourth class. The magnetic field is a small 
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Figure 4.9: Almost non-magnetic state m = 10, n = 62, B = 25. The numerical results are 
shown as before, except for the upper-left corner, where the theoretical streamlines are shown. 
The theory and numerics disagree mainly in the low-current regions. The diamagnetic loops 
close to the diagonals nearly cancel the paramagnetic loops in the triangular wedges between 
the diagonals. 
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Figure 4.10: Normal derivative \d^ nm /dy\ y= _i/2, obtained numerically (upper plot), and the 
theoretical surface of section function ^ m {x) (lower plot) for the state of Fig. 4.9. When m 
is large, ip m {x) extends considerably outside the domain [—5, 5] (indicated by the vertical 
lines). In order to construct the normal derivative, the tails have to be "folded" back with an 
additional phase factor. For small m both functions are proportional. 
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Figure 4.11: Numerical wavefunction and current for the state of maximum diamagnetism 
m = 14, n = 62, B = 25. The same representation as Fig. 4.5. 
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Figure 4.12: Theoretical wavefunction and current for the state in Fig. 4.11. The same repre- 
sentation as Fig. 4.5. 
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Figure 4.13: Numerical wavefunction and current for the diamagnetic state m = 14, n = 62, 
B = 31.4. The same representation as Fig. 4.5. 
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Figure 4.14: Streamlines for the sequence of states m = 6, 8, 10, 12, n — 62, B = 25, counter- 
clockwise from the lower right. 
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perturbation in this case and the states are weakly diamagnetic. The expansion near (1, 1) 
resonance will eventually break down and higher resonances become involved. 



4.2 Aharonov-Bohm flux line 

The theory laid out in the previous sections can be extended to non-uniform flux configurations. 
If the Lorentz force is neglected then the problem is reduced to finding the flux enclosed by a 
periodic orbit. Of course, if the flux distribution is not smooth enough, like the ABFL case, 
we have to worry about the diffraction. To estimate the effect of diffraction consider a finite 
size flux line of radius p and flux <f>. An orbit that goes through the flux region is subjected to 
the Lorentz force produced by the magnetic field B = <fi/irp 2 and gets deflected by an angle 

(4.20) 



kp 

where 4>q is the flux quantum. The diffraction can be neglected if 69 -C 1. Clearly, in the 
zero radius flux line the diffraction effects cannot be neglected. If we define the perturbation 
parameter e = (0/0o) /k, i.e. the ratio of the perturbed part of the action to the unperturbed 
action, then the condition is e/p <C 1. This means that the limit p — > cannot be taken before 
e — > if we want to neglect the diffraction. Therefore we have to consider a finite size flux line. 
In the numerical examples we take <p/4>o — 0.1 and p = 0.01, while k > 120. Instead of one 
flux line, four quarter-strength lines symmetrically located were used in order to reduce the 
amount of calculations. Within the approximation of our theory it gives essentially the same 
result as a single line. Note that for the ideal ABFL the matrix elements of the Hamiltonian 
are infinite if the basis functions do not vanish on the flux line. One can avoid this problem by 
making the substitution ^ — > T-l'/'/^ol^r w here r is the distance from the flux line |56|, |6l] . The 
new wavefunction ^ satisfies the equation H^> = E^> with the same energy as The most 
singular term of the non-hermitian Hamiltonian 7i is of order r _1 instead of r~ 2 . 

We consider the states near the (1, 1) resonance. A periodic orbit of this family may either 
enclose the flux line and have ±27r0/</>o added to its action, or pass by the flux line and have 
no extra action. Hence the one-dimensional wavefunction satisfies the Schrodinger equation 



(4.2) with a step-wise effective potential. Strictly speaking, the walls of the potential will have 



a finite width p, because when the orbit goes through the flux line there is an additional action 



119 



between and ±2ir<p/<po- Of course, the profile of the walls depends on the flux distribution 
within the line. If we are not going beyond the second order of the perturbation theory (Sec. 
|2.2| ), we may disregard the finite width of the walls, since the wavefunction does not change 
much. This way the results would be independent of p. The third order contains the derivative 
of the potential. In order to neglect it we must require ke 3 ^ 2 /p <C 1. Since e/p <C 1, we have 
the condition k^/e < 1, or ke -C 1. 

The width and position of the square well and barrier depend on the location of the flux. 
If the flux is located at x = 0, y = — \ + a, where < a < 1/2, the potential is 



V (x) = < 



~ 2 *&b xe[-a,a] 

+2tt££, xe[l-a,l + a] (4.21) 
0, x € [a, 1 — a] 



extended periodically by V (x + 2) = V (x) . Hence, assuming <p > 0, the potential consists of 
the well of width 2a and depth 2-k-^-^- in the interval [—5, 5] (counter-clockwise orbits) and 
the same size barrier in [i, |] (clockwise orbits). For sufficiently large the eigenfunctions 
are approximately ip m (x) = cos [(to + 1) irx/2a + to.7t/2] , |a;| < a, and zero elsewhere. This 
expression holds for sufficiently small m. The energy E m w — 2ir-^jk + v ^t"^ • 



The function ipo for a — 1/4 and n = 86 is shown in Fig. 4.15 (a). We include the 
exponentially small portion of the function for |a;| > 1/4. The upper curve is ipo (%) and the 
lower curve is ipo (—1 — a;) . The PSS wavefunction ip (x) — e lKX ip (x) should be compared to the 
normal derivative of the two-dimensional wavefunction |95 , /9n| at y = —1/2 [curve (b), Eq. 



(4.8)]. The oscillations appear because of interference of ip (x) and ip (—1 — x) , they scale as 
ipo (— 1 — x) . The other four curves show the normal derivative calculated numerically. Curve 
(c) is from full numerical diagonalization, (d) is from the diagonalization in the reduced basis 



that includes only the unperturbed states along the diagonal p=j + l,q=j~ I (Sec. 4.1.4 ). 
The latter result is very close to the theoretical. Curves (e) and (f) are the numerical results for 
the ideal single ABFL (p — 0) with n = 82 and n = 70. Although there are strong diffraction 
corrections, the overall shape is given well by the theory. The two-dimensional wavefunction 
is localized along the stable orbit. The strength of localization is characterized by ip$. Figure 



4.16 shows I 'J' (x,y)\ for an ideal ABFL, n = 58. Although the localization is not very strong, 
the wavefunction has little support in the center and the corners of the square. 
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Figure 4.15: State m = 0, n = 86 in the square with the flux line, a = |, <f> = O.l^o, P = 0.01. 
(a) Upper curve ^po(x), lower curve ^{—l — x). This function is superimposed on the other 
plots that show the normal derivative /dn\ y= _i/2- (b) is the theoretical result from Eq. 
(4.8), (c) is from the full numerical diagonalization, (d) is from the diagonalization in the 
reduced basis that consists of the states along the diagonal p^^ + Z, <Z = § — I (Sec. 4.1.4). 
For comparison, curves (e) and (f) show the numerical results for the ideal single flux line 
with n = 82 and n = 70, respectively. Although there are strong diffraction effects, the overall 
shape is given correctly by the theory. 
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Figure 4.16: Absolute value of the wavefunction, n = 58, m = 0, for an ideal single flux line, 
a = |, = O.l^o- Note that |^(a;,y)| is almost but not strictly symmetric under the 90° 
rotation due to the higher order effects. 
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If the flux is located at the center of the square (a = 1/2), the potential well extends to 
the whole side of the square and there is no spatial localization. For small m the counter- 
clockwise orbits arc preferred and there is a strong paramagnetic current. The upper part 



of Fig. 4.17 shows the streamlines for the state n = 82, m = 0, symmetric under the 
90° rotation. The lower part of the figure gives the current density along the line y = 0, 
x G [— 5,0] . Neglecting the A-term in Eq. (4.18) we have an approximate formula for it, 



j y (x, 0) cx — (cos7mx/2) 2 ip m (— x — — tp m (x — 5) . For low m the second term in 
the square brackets is small when x < 0, so j y < 0. It has equally spaced double zeros at 
x = (21 + 1) /n. The factor depending on i\) m has a zero at x — 0. 

It is interesting to compare the square with ABFL with the "step" billiard. This is a square 
with the step-wise boundary perturbation. For example, suppose the lower side is perturbed 
by £ (%) = — e, |x| < a, and £ (x) = 0, \ > \x\ > a. Then for the (0,1) resonance (bouncing 



ball states with the large y-momentum) the effective potential is proportional to Eq. ( 4.21 ) 
limited to x E [ — -i , -i] with the periodicity V (x + 1) = V (x) (the system has the time- reversal 
symmetry). When fee <C 1 and k^fl ~ 1, our theory predicts localization within |x| < a. If 
ke = 7r/2, i.e. 2e is a half- wavelength, the classical action for the (0,1) orbits within |x| < a 
and outside of this region differ by it. This is equivalent to the case 4>/4>o = 1/2- Since the 
action enters the phase of the T-operator, this difference becomes ambiguous, it can be made 
— 7r by adding the phase 2-7T for |x| > a. (This would be impossible in the case of continuous 
perturbation.) Thus the effective potential is undetermined when ke > 1. The preliminary 
numerical results J6l[ show that in this case there are states localized within |x| < a and the 
states localized outside of this region. 

4.3 Experimental suggestions 



Some of the experimental methods mentioned in Sec. 2.3.6 could conceivably be adapted to 
the billiard with a magnetic flux. The mesoscopic systems like the quantum corrals [p7| or the 
GaAs squares f3lf are directly related to our theoretical model, though accurately measuring 
the wavefunction may be a challenge. 

In the liquid surface wave experiments the effect of the flux can be modeled by the flow of 
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Figure 4.17: Ideal flux line at the center of the square, cf> = O.lfa. Current in the state n = 82, 
m = 0, symmetric under the 90° rotation. Upper figure: numerical streamlines for a quarter 
of the square. Lower figure: theoretical current density j y (x, 0), dashes indicate the numerical 
minima. 



124 



the medium |T(J . If V (r) is the velocity of the liquid then there is a correspondence 

eA(r) -fcV(r) 



(4.22) 



ch v g (k,r) 

where the group velocity v g <C |V| . (This analogy neglects the A 2 term in the Hamiltonian.) 
Thus, rotating the tank with a constant angular velocity is equivalent to applying the uniform 
field to a quantum billiard. The flux 

2tt-^ < — t/V-dr (4.23) 
<Po v g J 

in the case of homogeneous medium. A flux line would be analogous to a vortex formed by 
the water pouring through a small hole on the bottom of the tank. Reference reports the 
experiments on scattering of the surface waves on a vortex. 

In the microwave field experiments |73, pq, ^2) a ferrite strip of length 2a can be embedded 



in the wall y = —1/2 (Fig. 4.18 ). If one applies the static magnetic field to the ferrite, the 
wave will acquire an additional phase upon reflection from the strip. The phase is different 
for the forward and backward directions, i.e. the time reversal symmetry is broken. This 
system is analogous to the ABFL square billiard with the flux line located distance a from the 
bottom. Suppose the square cavity is in xy-plane, and the ferrite is under the magnetic field 
in z-direction. Then the permeability of the ferrite in the absence of losses is 



(4.24) 



Consider a plane wave with electric field ~E, i =Eie l ( kxX ~ kyy ^ (k y > 0) propagating in the cavity. 
Suppose this wave is incident on the ferrite layer of width S. Assuming the Dirichlet conditions 
on the metal wall to which this layer is attached, one can show that upon reflection the wave 
becomes E r = —Eze l ^ kj:X+kyy+ '^ where the phasef] 

k y (jj, 2 - x 2 ^j sin (k^S) 
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2 tan 



fcy/tty cos (ky8^ — k x xs\n (k{j&^j 



(4.25) 



Here k^ = yk 2 — k 2 where kf = fcy^uj| — >* 2 /m|| i s the wavenumber inside the ferrite. The 
time-reversal is achieved by reversing the sign of k x . Then <f> changes. In the quantum billiard 



This differs from the result reported in Rcf. |72| 
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the phase changes sign but keeps the magnitude when time is reversed. Here it is not the 
case because <fi includes the phase that the wave acquires by traveling through the bulk of the 
ferrite. In order to separate the effect of the magnetic field from the effect of the width we 



may remove a metal layer from the wall outside of the strip, i.e. for |x| > a (Fig. 4 . 1 8| ) . The 
width of the removed layer 6' is determined by the condition 

2k v S , = ^[<f>(k x ) + 4>(-k x )] (4.26) 

meaning that the phase the wave acquires outside of the ferrite strip should be equal to the 
average phase on the strip. Then the effective phase <p e g — <f> — 2k y S' is odd in k x . In the 
experiments the wavenumber in the ferrite kf was about ten times larger than k, so we expect 
5' > S to compensate for the phase. 



4.4 Conclusions 

The classically weak magnetic field acts as a perturbation in an integrable square billiard. 
In this system there are special states that are localized near the short stable periodic orbits 
that enclose a finite flux. The perturbation breaks the time-reversal symmetry and there is a 
preferred current direction. Not surprisingly, the states carry persistent currents. We observed 
a variety of probability and current distributions. The overall magnetic response ranges from 
the paramagnetic for the well-localized states to the diamagnetic. The short periodic orbits 
dominate the susceptibility. The localization also takes place for an off-center Aharonov-Bohm 
flux line. Although there is no direct classical effect in this case, the vector potential changes 
the action and enters the phase of the wavefunction. The diffraction effects are quite important 
for the flux line and limit the preciseness of the semiclassical results. 
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y =1/2 i 




removed layer 



Figure 4.18: Possible microwave experiment with ferrite. The ferrite layer of width 8 and 
length 2a is mounted on the wall. Part of the wall of width 5' is removed to compensate for 
the width of the ferrite. The dashed line indicates the boundary of the original square. 
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Chapter 5 

Quasiclassical Born-Oppenheimer 
approximations 



Up to now all our results followed from the Bogomolny equation. It gives the natural semiclas- 
sical formulation and reduces the dimensionality of the problem through the use of Poincare's 
surface of section (PSS). The perturbation theory allows one to make an asymptotic expansion 
in the phase of the wavefunction in a controllable way and relate the quantum states to the 
classical phase space. The method has certain shortcomings. To find the full two-dimensional 
wavefunction an additional trivial, but lengthy, step is required. The resulting function is 
inaccurate near the caustics and in the classically forbidden regions. 

There are cases, on the other hand, when the equivalent results can be obtained by the adi- 
abatic, or Born-Oppenheimer, approximation (BOA) p2| . In general, the BOA is appropriate 
when a subset of the system's coordinates varies in time slower than the remaining coordi- 
nates. This allows an approximate separation of variables in the partial differential equation 
describing the system. We will see in the examples below that such separation of variables is 



possible sometimes even if the notion of "fast" and "slow" is not well defined (Sec. 5.5). In 
any case, the conditions for the separation always follow from the differential equation itself, 
and we refer to all such cases as the BOA. 

The relationship between the BOA and the T-operator method is quite complex | p0| . Both 
methods are usually equivalent in the leading order. Unlike the T-operator, the BOA directly 
generates the two-dimensional wavefunction. If the PSS is chosen along the slow direction, the 
one-dimensional PSS wavefunction will be a part of the BOA result. The T-operator method 
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often produces an asymptotic form of the BOA solution. It is not always clear whether the 
BOA can be constructed for a given problem. For instance, we were unable to formulate the 
BOA in a simple form for a perturbed circle, apart from the states near the (1,2) resonance 



(Sec. 5.4). The T-operator, once constructed, treats all resonances uniformly. The purpose of 
this chapter is to illustrate the above remarks on several examples. 



5.1 Textbook example 



We remind the reader of the standard example where the BOA can be used In the 
molecules or solids the electrons' positions r e can be treated as "fast" compared to the "slow" 
ionic positions R^. This is based on the small electron- ion mass ratio m e /Mj. In the Schrodinger 
equation 



2m e e 2M, t 
we make the Born-Oppcnhcimcr ansatz 



V 2 + ^(r e ,R0 



*(r e ,Ri) =E^ (r e ,Ri) 



*(r e ,R() = $(r e |R 2 )^(R,). 



(5.1) 



(5.2) 



Here we assume that $ (r e |Rj) is the Nth electronic eigenstate for fixed ionic variables which 
solves 



2m, 



■Vl + V(r e ,Ri 



$(r e |R 2 ) = £/(R,)$(r e |R 4 



(5.3) 



where Ri is treated as a parameter. The eigenenergy U (R,) then acts as a potential for the 
slow variable: 



2M 



<V? + t/(Ri) 



if> (Ri) = EiP (Ri) 



(5.4) 



The adiabatic invariance provides that the electronic eigenstate label N does not change as 
Ri is slowly varied. The approximation made to the Schrodinger equation implies that 



ViV>Vi$ 



Mi 



V>V?$ 



Mi 



(5.5) 



In the case of the ground state of hydrogen molecule one can estimate |V,$| ~ |V e $| 



|$| /as where as is the Bohr radius. The ions will oscillate with frequency uj ~ e / -J ' Mia A B 
and amplitude SRi ~ h/Muj ~ (rrie/Mi) 1 ^ 4 as- Hence we estimate |V^| ~ \ip\/5Ri ~ 



\tp\ (Mi/me) 1 / 4 /olb- Therefore the first term on the left of Eq. (5.5) is the biggest and it 
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is (Afi/m e ) 3 ^ 4 times smaller than the r.h.s. This example shows that the BOA amounts to 
neglecting certain derivatives in the differential equation. The higher order corrections can 
also be written down. 

It was shown that the interaction between the fast quantum electronic and the slow classical 
ionic degrees of freedom may lead to chaotic behavior |L7j. 



5.2 Bouncing ball states 

In Sec. |3 - 2) we have shown how to derive the bouncing ball (BB) states in a Bunimovich stadium 



(Fig. 3.1C) using the T-operator. We obtained an almost square well effective potential which 
made the states localized within the straight part of the billiard. Equivalently, we could use 
the BOA H since there is a separation between the fast and slow motion. We have also 
mentioned a generalization of this problem: the stadium with the endcaps of slightly different 



radii (Fig. 3.11). In the latter case the effective potential is a square well with a sloped bottom, 
and the states are spatially shifted towards the wider part of the billiard. In both cases the 
wavefunction stays away from the semicircular regions, so the exact shape of the boundary 
there is of little consequence for the states with the low transverse quantum number m. Another 
interesting example is the 7r/3- rhombus billiard fus)] , which has two degenerate families of the 
BB states. 

In this section we apply the BOA method to the (0, 1) resonance states in the tilted unit 



square introduced in Sec. 3.1.4. The lower side of the square has a slope, y = ex, and the states 
of interest have a small x-momentum compared to the y-momentum. These are, of course, the 
BB modes, and, as we just explained, they are similar to the BB states in a tilted stadium 
(the effective potential is twice as large in the latter case, because two sides are tilted). 

We solve the Helmholtz equation (V 2 + k 2 ) * = with the BOA ansatz 'J' (x, y) = 



<& {y\x) -0 (x) . The fast equation is 



<9 2 $ 
dy 



— = -U{x)<P. (5.6) 



The solution 



.r> ( </|.r. X i^^fa{ m ±-JL) ,5.7i 
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vanishes at y = ex and y — 1 and is normalized. From this we find the function U (x) = 
Inn I (1 - ex)] w TT 2 n 2 + 2e7T 2 n 2 x. Then the slow function satisfies the equation 

- ip" + a 3 xip = £ m tfj (5.8) 

1 /3 

where £ m = k 2 — 7r 2 n 2 and a = (2e7r 2 n 2 ) . Note that a 3 x is in the leading order k 2 e times 
the effective potential 2x of the T-operator method if the lower side is chosen as the PSS. 
Hence tp (x) is the surface of section wavefunction given by Eq. (3.11) with E m = £ m /k 2 e. The 



energy k 2 and the two-dimensional wavefunction are also given correctly. 



Equation (5.8) has a solution 

tp (x) = Ai (ax — z m ) (5-9) 

for large a. Here z m is a root of Ai (— z), so that ip (0) = 0. If a — z m ^> 1, ip also effectively 
vanishes at x — 1. The transverse energy £ m = a 2 z m must be much smaller than the maximum 

1/3 

potential a 3 , i.e. a ~ {k 2 e) 3> z rn . The wavefunction ^ (x,y) is localized within x < z m /a. 
If this condition is not satisfied, Eq. ( |5.8| ) still can be solved by a linear combination of Ai 
and Bi functions. In the opposite case, k 2 e <C 1, the slope of the side can be neglected in the 
leading order. This is the situation when the non-resonant perturbation theory applies and 
there is no localization. The BOA is applicable if j*" 1 (d$/dx)\ |V>"V| < I* -1 (d 2 <5>/dy 2 ) \ . 
When z m /a -C 1, we estimate ~ a : an d the above condition gives e 4 -C k, which is 

always true. For larger m, i.e. when z m /a ~ (ra/k^/e) 2 ^ > 1, we have ~ m i so the 

condition becomes em <C k ~ n. This shows that the BOA is valid even if m ~ n. 



Figure 5.1 shows the cross-sections of the numerically obtained wavefunctions ^ nm ■ (a) 
^55,1(^ = 0.01,2/) [this is proportional to $ (y\x = 0.01)], (b) ^55,1 (a;, y = 0.99) compared 
with Ai (ax — z\) , and (c) ^55,2 (%, V = 0.99) compared with Ai (ax — 22) ■ Here e = 0.01 is 
less than the wavelength A = 0.036 which in turn is less than ^fe = 0.1. This results in 
sufficiently strong localization. The two-dimensional representation of | ^55,2 (x,y)\ 2 is shown 



in Fig. 5.2 



One concludes that in the case of the (0, 1) resonance in the tilted square the BOA is 
straightforward and readily provides the two-dimensional wavefunction, while the T-operator 
method would require more work. The BOA method is not, however, directly generalizable to 
the higher resonances (Sec. [3.1.4 ), although the next section suggests a possible approach. 
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Figure 5.1: Cross-sections of the numerical wavefunctions \&„ m for the (0,1) resonance in 
the tilted square, e = 0.01. (a) '£55,1(2; = 0.01, y); (b) *55,i(a;,y = 0.99) compared with 
Ai(ax — zi); (c) & 55,2(2, y = 0.99) compared with A\(ax — z 2 ). The magnitude of the states 
has been normalized. 
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Figure 5.2: Density plot of 1^55,2(2;, y)\ 2 , the state of Fig. 5.1 (c). The lower part of the billiard 
(0 < y < 0.25) is shown and the graph is expanded in the y direction by a factor of three in 
order to display more details of the wavefunction. The dashed line is at y = 0. 
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5.3 Channeling approximation 

The necessity to find the coordinates in which the fast and slow motion separate may prevent 
us from treating generic resonances by the BOA. In some cases the method of images resolves 
this problem. Consider, for example, a rectangular billiard perturbed by a potential or a 
magnetic flux. The copies of the rectangle obtained by reflection about its sides will cover the 
whole two-dimensional plane, forming a lattice. (The magnetic flux has opposite signs in the 
neighboring rectangles.) Suppose that the perturbation is classically weak, that is the classical 
orbit does not deviate much from the straight line after one passage across the billiard. (For 
the estimate in uniform magnetic field see Sec. [4.1. l[ ) An orbit in the original rectangle extends 



along a straight line with minor deviations in the extended scheme (cf. Fig. 4.1). Now the fast 
and slow directions separate and the BOA can be used. The BOA in this case is related to 
the channeling method for the energetic particles traveling through the crystal lattice j5^| . 

For definitiveness, consider the uniform magnetic field B as a perturbation. We start with 
the (1,1) orbits. We rotate the coordinate system to make one variable £ = (ax + by) /d 
run along the fast (1, 1) direction and the other variable r\ = (— bx + ay) jd along the slow 



perpendicular direction. Here a x b are the dimensions of the rectangle and d = V 'a 2 + b 2 . The 
Schrodinger equation in these coordinates in the units h = c = e — 2m = 1 is 

{ [-id € ~ A,: (£, r,)] 2 + {-id v - A, (£, r?)] 2 } * (£, v ) = E * (£, r,) . (5.10) 

In the channeling approximation we assume that the particle moving fast in the £ direction 
"sees" only the averaged vector potential. We may neglect A 2 comparing to the linear kA 
terms (since e = B/k <C 1) and substitute the average vector potential A^ (rj) for A^ (£,??) , 
where 

As(r ] )=L- 1 d^A e (C,v) (5.11) 
Jo 

and L(: = 2d, the period in £. Note that A^L^ = § A ■ dl is the flux enclosed by the periodic 
orbit in the original billiard. Similarly, A n L^ = § dl x A = in the gauge where divA = 0. 



Now Eq. ( 5.10 ) becomes separable. 

The same result can be obtained on a more formal basis. Namely, we make an ansatz 
^ (£j 7 /) = ^ iClv) "0 {v) an( i consider the vector potential as a small perturbation in the fast 
equation for Then, according to the elementary perturbation theory, we can approximate 
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<I>(£|?7) = e lk ^ and find the fast eigenenergy U (77) after multiplying the fast equation by <&* 
and integrating over £. Again neglecting A 2 we find that 

U (77) = k 2 - 2kA^ (77) . (5.12) 

The A^c^ term automatically disappears in the BOA approach. The slow equation is 

- 7//' - 2kAf. (77) t> = (E - k 2 ) il). (5.13) 

Now we can rewrite the results in the x,y coordinates. Suppose the origin is chosen at the 
center. Define the function ip (x — ay/b — a/2) = ip (— drj/b — a/2) = ip (77) . This definition 



relates tp to the function in Eq. (4.7) when a = b = 1. The above equation turns into Eq. (kL3), 



fa (x) + [E m - V (x)} i> m (x) = 0, (5.14) 

with V(x) = -2$gkAt(-!f - . It has the symmetries V (x) = V (-x) and 
V (—a/2 + x) = —V {—a/2 — x) and has a period 2a. For uniform field 

= < °^^> (5.15) 

with A ? (77 + 2ab/d) = A\ (77) , and V (x) = - (Bkab 3 /d 3 ) \-2 (§) 2 for |x| < a/2 repeated 
antiperiodically outside this region. An eigenfunction in this potential satisfies the Bloch 
condition 

4> m (x + 2a) = e i0 i) m (x) (5.16) 

where (5 is to be determined. 

The two-dimensional wavefunction has the form 

* (x, 77) = e lk( - ax+b ^' d 4> m (x - ay/b - a/2) . (5.17) 

This is one of the four degenerate solutions. Another solution is obtained from this af- 
ter the rotation by tt, ^2 (x, y) = ^0 ( — x, —y) ■ The other pair comes from shifting to the 
neighboring square, where the magnetic field has the opposite direction, and considering the 
(—1, 1) orbit. These solutions are \&i (x, y) = exp [ik (—ax + by) /d] ip m (— x — ay/b + a/2) and 
^3 (x, y) = $1 (—x, —y) . The wavefunction in the original rectangle is constructed from these 
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four solutions under the condition that it vanishes on the boundary. The state can be made 
even or odd under (x, y) — > (— x, —y) . With r — ±1 it has the form 

y nm =A(y + ry 3 ) + B(y 2 + ry±). (5.18) 

Requiring * nm (±| ,y) = 0, we find that A = -Be~ ix " and e i/3 = e~ l2x " where Xa = ka 2 /d. 
And from $ nm (x, ±|) = we get A = —rBe lXb and e J/3 = e l2xb where \b = kb 2 /d. Therefore 
Xa + Xb — where n is even for r = 1 and odd for r = — 1. This determines /c = nn/d, and 
the energy can be expressed as 

E mn = (irn/d) 2 + d 2 E m /b 2 . (5.19) 

The phase shift j3 — 2Trn(b/d) 2 — 2imb- If we assume \/3\ < ir then n is an integer part of 
(nb 2 jd 2 — 1/2) . For a square f3 = for even n or it for odd n. This is the case considered in 
Ch. [|. It is clear that the BOA and the T-operator theory agree at least to the order ^ft in 
the wavefunction and e in energy. 

We were able to apply the BOA to the rectangular billiard by extending it to a bigger 
system and thus "unfolding" the trajectories. The method can be generalized to include the 
boundary perturbation. In this case we again extend the unperturbed rectangle by reflection 
to the whole plane and assume that the classical orbits arc almost straight lines in the plane 
(actually, they will slightly bend on the boundary crossings). Then U (77) — k 2 — 2k 2 SL (rj) /Lg, 
where SL (rj) is the additional path due to perturbed boundary for an unperturbed orbit of 



period L^. In light of Eq. (5.10) we might say that the change in length is compensated by 
a change of momentum 5k = —k5L/L that preserves the action. Notice also that the time it 
takes the particle to travel distance L with the modified momentum is equal to the time it 
travels distance L + SL with the old momentum. If the extension of a system to the whole 
plane in the original coordinates is not possible (like for a circle billiard) , it can be done in the 
action-angle variables. 



5.4 BOA in the asymptotic region 

Now we return to the low angular momentum resonance in a perturbed unit circle (Sec. [2.2.2 ). 
In this example the separation of fast and slow coordinates holds only over a part of the 
billiard. The BOA wavefunction in this region is in agreement with the T-operator solution. 
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The Hclmholtz equation in cylindrical coordinates has the form 



9 8 , 22 
r— r— + k z r z + 
Or Or 



d 2 



* = o. 



(5.20) 



For the low angular momentum I states the angular coordinate is slow and the radial coordinate 
is fast everywhere except the region near the center r < l/k <~ \fk. Outside of this area we can 
make the BOA ansatz (r, 9) = $ (r\6) tp (9) . The fast equation 



r— r— $ + [fcV - I 2 (0)1 $ = 
Or Or 1 J 



(5.21) 



is the Bessel equation with the variable order I (6) . Since it is invalid for r — > the solution 
may include both Bessel and Neumann functions. It is convenient to invoke the asymptotic 
expansion for large kr and / but small I / kr EJ] 



$(r|6») 



1 



kr 



kr - 



2kr 



+ a(6) 



(5.22) 



where a (9) is the phase that mixes Bessel and Neumann functions. Function $ must vanish 
at the boundary r (9) = 1 + eAR (9) which requires 



;2 iq\ _ 1 

fc + /ceAi? (6») + — K -l — I +a(9) = Tr[n-^- 



2k 



(5.23) 



Both functions I (9) and a (9) are unknown, and we need some additional information to 
determine them. Consider the slow equation 

V>" + I 2 (9) i> = 0. (5.24) 

The radial wavefunction ip (9) is expected to coincide with the PSS wavefunction that solves 



Bogomolny's equation. Equation (2.10) implies that 

I 2 (9) = k 2 e [E m - V (9)} 



(5.25) 



where V (9) = AR (9) + AR (6 - it) . Since V (9) has period 7T, Eq. Q5.24| ) is solved with the 
boundary condition iJj (9 — it) = (—l) m Tp(9), and the eigenvalues £7 m are determined. Now 



Eq. ( 5.23 ) enables us to find 

a (9) - 

and the quantization condition 



-fee [AR {9) -AR(6- tt)] + a 



a = tt ( 71 - - 



(5.26) 



(5.27) 
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The constant ao is determined from the special case AR = 0. Then I = m and $ (r\9) oc 



Jm (kr) , which implies ao = — ^ (mmod 2) — j. Figures 2.5, 2.6 show the states with m — 2, 3, 
respectively, in the short stadium. 

In the considered example the BOA is not self-contained — we could not complete the 
solution without knowing the behavior of the wavefunction at the origin. Having obtained the 
radial wavefunction by other means (e.g. from Bogomolny's equation), we reconstructed the 
two-dimensional wavefunction [cf. Eq. ( 2.6l| )] and the quantization condition [cf. Eq. (2.11)] of 
the perturbation theory. It is remarkable that the fast wavefunction <f> (r\9) becomes asymp- 
totically the Bessel function Jug) (kr) only in an unperturbed circle, in all other cases it is 
mixed with the Neumann function N^g^ (kr) . 



In this section the phase of the wavefunction ( |5.22j ) and the quantization condition ( 5.25 ) 
contain a k~ l order corrections that were absent in the T-operator theory. This corrections 
are normally small in the semiclassical regime although they may begin to play role when 
the angular momentum I ~ ky/e < 1. There are several sources of the A: -1 corrections in 
the T-operator method. First, the terms can be added to the phase of the T-operator 
itself. For example, in the boundary integral method these terms originate in the derivative of 



Hankel's function [cf. Eq. ( 1.13 )1- Second, the Bogomolny equation ip = Tip should be solved 
to higher precision, namely, one has to retain the next order terms in the expansion of the 



unperturbed action and the prefactor near the stationary point [cf. Eq. (2.6)]. Similarly, the 



two-dimensional wavefunction can be found from Eq. (1.10) with the Hankel function as the 



kernel. Note that the more general kernel (2.54) is valid only in the leading order in fc _1 . 

The low angular momentum resonance in the perturbed sphere of Sec. |2.5.4 can be described 
by a similar technique. We start with a three-dimensional Helmholtz equation 

|_ r 2|_ + fc 2 r 2_^ = () (5 28) 

or or J 
where I 2 = — ^-q4e sin 0-^ — Ars-is- The equation in the fast radial direction 

— r 2 — \k 2 r 2 -1(1 + 1)] $ = (5.29) 
or or L J 

is an equation for a spherical Bessel function of order I = 1(9, ip) fil]] . Again, since the BOA 
is not valid near the center, the solution will include both Bessel and Neumann spherical 
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functions. Asymptotically, 



<& (r\9, ip) ~ -— sin 
kr 



where a (0, 95) is the phase to be determined. The slow equation 

- Pip + I {I + 1) t/j = 



(5.30) 



(5.31) 



implies that 

Z(Z + l) = fc 2 e[£; iro -F(0,^)], (5.32) 

where 1/ (0, 99) = AR (9, ip) + AR (tt — 9, ip + ir) . The eigenfunctions ipi m [the quantum num- 
ber I should not be confused with the function I (9, ip)] will be even or odd under inversion, 
ipi m (9, ip) = (—1)' ipim (n — 9, ip + tt) , where we assume the state tpi m has the same symmetry 
as the unperturbed state Yi m (9, <p) . As before, the Dirichlet condition for $ gives 



a{9,ip) = --he [AR{9,ip) - AR (tt - 6, ip + tt)] - - (/mod 2) 



and the quantization for k 



1 TT 

k + -keEi m — TTfi + — (I mod 2) . 



(5.33) 



(5.34) 



Clearly, the BOA is consistent with the T-operator approach and even improves on the angular 
differential equation. 



5.5 Whispering gallery modes 



The whispering gallery modes discussed in Sec. 3.3 can also be derived with the BOA. The 
adiabaticity comes from the slow variation of curvature of the boundary. This approach is 
related to the parabolic equation and the etalon methods Q in the sense that in the former 
some derivatives in the partial differential equation are neglected and in the latter an ansatz 
involving the Bessel function with variable order is made. All these methods give the same 
leading order results. 

The billiard boundary is locally a circle of radius R (s) where the variable s runs along the 
perimeter. Therefore the whispering gallery wavefunction is locally a wavefunction for a circle. 

* (p, s) = $ (p\s) V (s) = a (s) J l(s) (At.) V («) , (5.35) 
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where p <§C i? is the distance from the boundary, r s — R(s) — p is the radius measured from 
the local center of curvature and a (s) is the normalization for The slowly changing angular 
momentum I (s) does not have to be integer and will be determined by the boundary conditions. 



The "slow" wavefunction ip (s) satisfies the equation similar to Eq. (5.24) and is equal to 

expi J s l(s')ds'/R{s') 



(5.36) 



where ds/R plays a role of dd. Since I ~ kR, the function ip (s) is not really slow, varying 
approximately as e lks . However, $ (p\s) depends only on the slow varying functions of s, like 
I (s) and R (s) , which justifies the BOA ansatz. Since the order of the Bessel function is large 
and close to its argument, it can be approximated by an Airy function. If we define a function 
f(s) byl = kR(l-f) then E| 



J l(s) {k[R(s)-p]} = J, 



1/3 f kRf'-kp 



kR 



1/3 



Ai 



,l/3 kR f - k P 



(kR) 



1/3 



(5.37) 



(I can be replaced by kR in the leading order in /' and p/R.) The Airy function must vanish 
at p = which makes 



kR 2V3 (fci?) 2/3 

where z n is the nth root of Ai (— z) . The full wavefunction can now be written down as 

' 2 l/3 fc 2/3 



(5.38) 



* (p, s) = (kR)- 1/6 Ai i^-^rj^P - z ") ^ ( s ) 



(5.39) 



where the normalization factor has been added. Note that the normalization ensures the same 
total current through any section s = const. 

The T-operator theory produces the asymptotic form of this solution, Eq. (3.41).^ There 
the function / (s) was defined via Eq. ( |3.33 ) as a part of the phase of the PSS wavefunction. 



Clearly, this is consistent with the current definition. The explicit expression for /', Eq. (B.34), 
compares with Eq. ( 5.38j ) if z n is used in its approximate form [|7r (n — ^)] 2 ^ 3 for large n. The 
BOA solution works both in the classically allowed and classically forbidden regions including 

1 Although denoted by the same symbol, the PSS wavefunction and the slow wavefunction differ by a 

factor R- 1 ^ 3 (s). 
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the caustic. It does not require the large number of wavelengths in the radial direction, which 
would be the condition for the asymptotic expansion of the Airy function. 

It is possible to derive the Airy function solution directly, without referring to the Bessel 
function. The Helmholtz equation in the (p, s) coordinates for small p/R is 
d 2 1 d 1 d 2 P^ 1 )' d 



* = 0. 



(5.40) 



dp 2 R(l-p/R)dp (l-p/R) 2 ds 2 {1-p/Rfds 
The second term can be neglected as compared to the first term since 

[according to Eq. ( |5.39| )1 and kR 3> 1. The third term can be expanded in p/R. The fourth 
term can be neglected compared with ^-§^z since ~ k and kR R' . The simplified 

equation 



1 + 2J R. 



dp 2 



k' 



ds 2 



* = o, 



(5.41) 



when solved by the BOA, yields Eq. ( 5.39| ). We can also estimate the limits of applicability 
of the BOA from this differential equation. When we make the BOA ansatz, the largest term 
that we drop is 



where we take the typical p ~ z n R 1 l z k~ 2 l z . We should compare it with 



<9 2 $ , 




k^ Ai" 






i? 2 /3 Ai 



Estimating |Ai"/Ai'| ~ n/z n we obtain the condition 

(fci?) 1/3 » \R!\z 2 Jn~ iR'ln 1 / 3 



(5.42) 



(5.43) 



(5.44) 



which is the same as the requirement (3.38) for the simple T-operator theory to work. In Sec. 



3.4 we consider the case when this condition is not satisfied. 



The classical caustic is given by the turning point of the radial equation, or a point where 
the argument of the Airy function vanishes, p(s) = ^^{^ ■ If the energy k 2 is fixed, the 
possible caustics are quantized by z n . In the classical picture the product e (s) R 1 ^ 3 (s) = I, 
where e (s) is an angle that the classical trajectory makes with the boundary, is an adiabatic 



invariant in the limit e — > (see Sec. 3.4). Semiclassically, I — 2 1 / 3 y / z^/fc 1 / 3 and, therefore, 
the caustics can be labeled solely by this parameter. This is not surprising since, apart from 
the energy which does not change the geometry of the orbits, this adiabatic invariant is the 
only (approximate) integral of motion. 
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5.6 Conclusions 



The Born-Oppcnhcimcr approximation can be used alongside the T-operator and other tech- 
niques to describe the perturbed integrable systems. When the separation of the fast and 
slow variables is possible it provides a convenient way to find a two-dimensional approximate 
wavefunction and energy levels. The BOA solution is usually valid near caustics and in the 
classically forbidden (shadow) regions where the semiclassical approximation often fails or 
needs modification. Whether the natural coordinates, in which the separation of fast and slow 
motion is possible, exist, depends on the geometry of the system or a family of orbits. In prin- 
ciple, the separation should always be possible in the action-angle variables, but then again 
the shadow regions are off limits. The adiabatic approximation may also be possible when 
there is a slowly changing parameter instead of the slowly changing position of the particle, as 
in the whispering gallery case. 
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Chapter 6 



Trace formulas 



The trace formulas are used to describe the density of states in terms of the classical periodic 
orbits of the system. The Gutzwiller trace formula jlO) applies to the hard chaotic systems 
and the Berry- Tabor trace formula (l|, [39| works for the integrable systems. In the former 
case the periodic orbits are unstable and isolated, while in the latter case they form families. 
The intermediate situations include the perturbed integrable systems J59fl and mixed systems. 

The T-operator formalism (Sec. |1.4| ) adequately describes a system on a semiclassical level. 
In particular, one should be able to derive the trace formulas directly from the T-operator. We 



show in Sec. 6.1 how this can be done. In the following sections we apply the general formalism 
to the Gutzwiller and Berry- Tabor cases. Then we consider a perturbed integrable system, 
which interpolates between the two extremes. The interpolation formula can be parametrized. 
At least four parameters are needed to describe the perturbed family of periodic orbits correctly 
|f79|], not three, as Ref. |76| claims. We illustrate these points with the example of the coupled 



quartic oscillators in Sec. 6.5 



6.1 General derivation of the trace formula 

The energy spectrum of a system described by the Bogomolny operator T(E) is given by the 



zeros of the Fredholm determinant [cf. Eq. (1.32)] 



D{E) = det [1 - T(E)] . (6.1) 
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The oscillatory part of the density of states can be expressed as a logarithmic derivative 



d osc (E) = d{E) - d{E) = — Im 

7T 



dlnD(E + ie) 



dE 



(6.2) 



where d(E) = '}2 a 8(E — E a ) and d(E) is the smoothed (Weyl) density of states. To justify 
this result |6(| we write D{E) = J[ n [l — e* 9 ™^)]. Here e 10 "^ are the eigenvalues of the 
T-operator which is approximately unitary. Then Eq. (|6.2|) yields 



dose 



-5E<( 1+1 "«'t)-E 



in 



(6.3) 



With 9' n approximately independent of n the first sum gives the smoothed part — d = —-^S' 
where N is the size of the T-matrix. The second sum is d(E). 

Using the relationship lndet(l — T) = Trln(l — T) in Eq. ( |6.2| ) and expanding the logarithm 
in powers of T, we find 

dose (E) = -Im y~] - dl ~ v ' E ' 



dE 



(6.4) 



where 



r n {E) = TrT n (E) = J dq x ■ ■ ■ dq n T \q x , q 2 ) ■ ■■T{q n ,q 1 ) 



(6.5) 



The composition property for the T-operator (Sec. 1.5) gives for a two-dimensional system 

S>S r {q,<i;E) 



\ 1/2 




J exp 

/ q=q' 


-S p (q,q;E) 



(6.6) 



Here p denotes a "closed" orbit that leaves the Poincare surface of section (PSS) from point q 
and arrives at the same point on its nth crossing of the PSS, S p (q, q; E) is the reduced action 
for this orbit, and the Maslov index was omitted. In effect, we calculated n — 1 out of n 
integrals in Eq. ( |6.5[ ) by the stationary phase (5$). Whether the last integral integral can be 
done by the as well, depends on the S p (q, q; E). 



Since the exact density of states is a collection of S- functions, the formal series (6T) diverges 
when E is on the spectrum. In practice, one includes only the first few terms that account 
for the short periodic orbits. This produces the density of states smoothed over some energy 
scale (7i/time). Such averaging can be experimentally relevant, say, due to the finite resolution 
of the spectrum because of non-zero temperature. For example, the susceptibility of a square 
in magnetic field is determined mostly by the shortest periodic orbit that encloses flux (Sec. 



4.1.3). Thus, only one term in the trace formula is needed in this case pj. When one is 
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interested in the energy level correlations on the scale of the average spacing, the orbits with 
periods up to the Heisenberg time hd should be included. 



6.2 Gutzwiller trace formula 

In the Gutzwiller case the periodic orbits are isolated, which means that S p (q,q; E) has the 
well defined stationary points. [dS p (q 7 q')/dq+dS p (q,q')/dq' = p(q)—p'(q') — at q = q' = q* } 
a S<& point of the integral ( |6 . 6[ ) . ] Let q, q 1 be in the vicinity of q* and expand in Sq = q — q* 
and Sq' — q' — q* 

S p (q, q') ~ S p (E) + p*(Sq - Sq') + ^V^Sq 2 + V 12 SqSq' + ^V 22 Sq' 2 (6.7) 

where p* is the momentum of the periodic orbit. The matrix of second derivatives V depends 
on q* and the energy. Suppose the orbit p of "length" n (i.e. returning to the surface of section 
n times) consists of r repetitions of a primitive periodic orbit of length s, with n = rs. Then 
the contribution of orbit p to the integral ( |6.6] ) is 

V 12 



E 



l/2\ 




J exp 


~{S P{ E) 



(6.8) 



\ l 

The summation is over the s stationary points, where the orbit crosses the PSS (each of them 
can be chosen as a starting/ending point for the periodic orbit). 

The prefactor can be expressed in terms of the monodromy matrix M p of orbit p, which 
relates the final momentum Sp — p — p* = VnSq + V\ 2 Sq' and position Sq to the initial 
Sp' = -Vi 2 Sq - V 22 Sq' and Sq' . With 





V^-VuVaa 


Via 


Via 


1 


V 22 


Via 


Via 



M p = Vl2 Vl2 (6.9) 



the prefactor can be written as s\ det(M p — 1)| -1 / 2 . The sum over the crossing points appears 
as a factor. [To see that det(M p — 1) is independent of the choice of the initial point q* , express 
M p = M 12 M 2 3 • • • M n i, where Ma+t is the monodromy matrix between the two consecutive 
crossings. Then det(M p - 1) = det(M 23 ■ ■ ■ M nl M 12 - 1).] 

For a repeated orbit, S p — rS s and M p — (M s ) r . Only the rapidly varying phase needs 



to be differentiated when evaluating the derivative in Eq. (6.4). With dS s {E)/dE = T s , the 
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period of the primitive orbit, the density of states 

T./h 



dose ^ 



|det(Mj- 1)| 1/2 
where we restore the Maslov index. 



cos < r 



S S (E) TT 

h 2 Vs 



(6.10) 



6.3 Berry- Tabor formula 

In the integrable systems the periodic orbits form continuous families, and therefore the 5$ 
cannot be applied. On the other hand, since the action is constant within a family, the integral 



(6.6) can be easily done in the action-angle variables. For a two-dimensional system we define 
the variables /, 9, J, 0, and choose the surface of section = 0. The action after the first 



return S(9 — 9' , E) depends only on the difference of the angles (cf. Sec. 1.7). For an orbit 
p of length n the action S P (A9) — nS(A6/n). If p is a periodic orbit with a frequency ratio 
uji/ujj — m/n then A8 = 9 — 9' = lixm. The derivative S" = — (2Trng'^)~ 1 > where J = gE{I) 



(Sec. 1.7). Using these results in Eqs. ( |6.6| ) and (6.4) we arrive to the Berry- Tabor trace 
formula @ [uj 



7T TT 

-v, 



(6.11) 



Note that a chaotic system may have a family of non- isolated periodic orbits. Among the 
examples are the bouncing ball orbits in the stadium |7^, W\ Q or Sinai billiard || (7(1 . 



6.4 Perturbed Berry- Tabor formula 



In a perturbed integrable system the families of periodic orbits are broken with only a few 



isolated orbits remaining (Sec. [2.l| ). However, the action S p (q, q'; E) in the integral (6.6) varies 
too slowly for the S<& to be applied. It may be convenient to express this one-dimensional 
integral in the action-angle variables. For an orbit of length n 



S p {0, 9') = nS{A9/n) + eW(6, 9') 



(6.12) 



where eW(9,9') is the perturbed part of the action that can be calculated by the standard 
technique. We may neglect the order e terms in the prefactor and take 9' = 9 — 2-kvh for a 
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broken family with the winding number m/n. Then the integral that remains is 



I W = — I (Wrxi) 



(6.13) 



where W{9) — W(9, 9 — 27rm).[|] The perturbed result for the d OS c is just the Berry- Tabor 
formula with the substitution cos0 — » Re[iw exp («</>)] where </> is the argument of cosine in 



Eq. (6.11) ]59[ . One might say that the perturbed trace formula interpolates between the 
Berry- Tabor case, when e — and Iw — 1, and the Gutzwiller case, when eW/h is large and 
I\v can be done by the S<&. 

When a family of periodic orbits is broken by the perturbation, at least one stable and 
one unstable periodic orbits remain, which means that W{9) has at least one minimum and 
one maximum. Assume, for example, that W(9) has a single maximum and minimum at 
9 = 0, 7r, respectively. Often the perturbed action W{9) can be approximated by a simple 
parametric expression that retains the essential information about its behavior. In the first 
attempts of this sort the function was approximated by the first terms of its Fourier series [^9| , 
e.g., W(9) ss wq + w\ cos{9 — 9 ). Then I\y exp(iewo/h)Jo(ewi/h) where J is the Bessel 
function. However, this expression is, in general, incorrect in the limit of large eW /h. Indeed, 



the integral (3.13) can be done by the S$ and, therefore, depends on the values of W and W" 
at the extrema. This means that Iw must depend on four parameters. 

Instead of increasing the number of terms in the Fourier expansion, it is more efficient to 
parametrize |?6| W(9) — Wq + W^cos^^)] where £(0) = and £(-7r) = 7r. The unknown 
function 0(£) can be approximated by the first two terms of its Fourier expansion 9 = £ — 



Asin£ - Bsin2£. With d9/dt; = 1 - A cos £ - 2Bcos2£ the integral (6.13) becomes 



Iw = ^r / d£(l- A cos £-25 cos 2£) exp 
2tt J 



— W\ cos£ 



(6.14) 



where Iq = exp[i(e/ft)Wo]. The integrals may be expressed in terms of Bessel functions 
J m (eWi/H), giving 

Iw = I (J -iAJ 1 +2BJ 2 ). (6.15) 
The four parameters Wo, W\, A, B, can be related to W and W" evaluated at the extrema 



1 Note that W(8) = rV q (9) for a periodic orbit (rp,rq), where V q (6) is just the effective potential of the 



T-operator perturbation theory (cf. Sec. 2.5.3). p, q are relatively prime. 
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by 

W(0) = Wo + Wi, W{tt) = W - Wi, 

Wi - Wi 
W (0) = — W"M = — (6 16) 

[ ' (l-A-25) 2 ' { ' {l + A-2Bf 1 ' 

Thus, the above interpolation formula is correct in the limit of large eW/H, as well as for 
e — > 0. One can show that det(M p — 1) is proportional to eW". 

The authors of Ref. |T6) (UGT) take B = 0. They compensate the lack of another parameter 



by letting the function g E [Eq. (6.11)] depend on e. Although formally this method seems to 
be correct, it is physically misleading. The function g E describes the constant energy surface 
in the action space (I, J). When the system is perturbed, the topology of the invariant tori 
changes and the old actions are no longer the integrals of motion. Therefore g" E does not 
have a clear physical interpretation for finite e. Secondly, UGT argue that "the independent 
evaluation of g E can be rather laborious and time consuming." However, in order to explicitly 
evaluate this function in the limit e — > (which must be equal to the standard g E ) UGT 
would need to know the parameter W\ (in our notation) and the values of det(M p — 1) for the 
periodic orbits in this limit, which in effect requires the solution of the unperturbed problem. 
With this solution on hand it should not be hard to find the standard g' E . 



6.5 Example: coupled quartic oscillators 

We study the correction to the trace formula for the coupled quartic oscillators, which are 
defined by the Hamiltonian 

H=^(pl+p 2 y ) + ax i + by 4 + ex 2 y 2 . (6.17) 

When e = we may introduce the action-angle variables (I, 9) in the x direction and (J, O) in 
y direction. For a given energy E and winding number a = Wj /wj 

1/4 



1 ( E/a 

a 
b 



V2 U + fa- 4 



2 r ( E/a 



sd(/c0), 

i=r^(jf^) 3/4 (6 - 18) 

and similar for and J. Here k = 2K/tt, where K = K(m = 1/2) is the complete elliptic 
integral of the first kind, and sd is one of the Jacobi elliptic functions El) . The constant energy 
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surface is described by 



9e{I) = J 



4/3 



1/3 



j4/3 



3/4 



The action after the first return to the PSS O = is 



S(6-6',E) = ^kE 3 / 4 



1/4 



(6.19) 



(6.20) 



a b 

The perturbed part of the action 8S = — J SHdt in the leading order in e. The integral of the 
perturbed Hamiltonian can be evaluated along the unperturbed orbit. For a family of periodic 
orbits of winding number a = m/n we find 

kE 3 / 4 a 



W{6) 



'%(ab) 1 / 4 [a + 6a 4 ] 3 / 4 



d6"s<r[n{6" + d)]sd 2 [K6"/a\ 



(6.21) 



2-7rm 



Note that the shape of W{9) is independent of a, b, and E, and each is taken to be unity. 



Figure 6.1 depicts W{9) for the (1, 1) family, a = 1. It is periodic with period tt due to 
symmetries. The interpolation formulas of the previous section can be adjusted, when the 
period of W(9) is 2ir/r, by changing 9 i— ► r9. From the numerical values of W and W" at 
the extrema we find the interpolation parameters A = —1.49 x 10~ 2 , B = —1.04 x 10 -4 
(independent of a, 6, E) and Wq = —0.46, W\ = 0.27. The real part of Iw as a function 
of e/h is shown in Fig. |3.2] and the imaginary part appears in Fig. |6.3|. The deviation of the 



interpolation formula from the exact integral is given in the inset of Fig. 3.2. They agree to 
within a part in 10~ J . 

Note that the parameter B is relatively small. This may have worked in favor of UGT 
who assumed B = 0. In general, however, this assumption implies a relation between the 
curvatures at minimum and maximum of W(9), namely, 



\/w 1 /w^ n + Vw 1 /\w^ x \ = 



(6.22) 



(in this case r = 2). We find numerically that the l.h.s. of this equation is equal to 1.00021 
for the (1, 1) resonance. 



6.6 Conclusions 

The oscillatory part of the density of states can be expressed in terms of traces of the powers of 
T-operator. The powers can be evaluated by the stationary phase. For their traces the 5$ can 
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Figure 6.1: W(9) for the (1,1) family of orbits (a = 1) of the quartic oscillators with x 2 y 2 
coupling. The parameters a = b = E = 1. 
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Figure 6.2: Real part of Iw as a function of e/h. The difference between this function computed 
by numerical quadrature and from the interpolation formula, |AJw|, is shown in the inset. 
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Figure 6.3: Imaginary part of Iw as a function of e/h. 
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be used only when the periodic orbits are isolated and the action of the closed orbits in their 
neighborhoods divided by Planck's constant changes by order of one on the scale shorter than 
the distance between periodic orbits (Gutzwiller case). For an (almost) integrable system the 
whole family of (almost) periodic orbits contributes to the trace (Berry- Tabor case). When 
the perturbed action becomes large compared to the Planck constant, the surviving periodic 
orbits can be considered well isolated and the Gutzwiller formula is recovered. 
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Chapter 7 



Summary 



In this work we formulated the semiclassical theory for perturbed integrable systems, applied 
this theory in a number of cases, and compared it with other methods. The use of the Poincare 
surface of section (PSS) allows to reduce the spatial dimensions by one. The Bogomolny T- 
opcrator connects the classical and semiclassical pictures on the surface of section. It is a kernel 
of an integral equation Tip = that determines the PSS wavefunction. It has a solution when 
the energy is on the spectrum, i.e. when det(l — T) = 0. The T-operator is unitary in the 
stationary phase approximation. When the boundary of a billiard is chosen as a surface of 
section, Bogomolny's equation follows from the boundary integral method in the semiclassical 
approximation. The T-operator contains the Maslov phases that come from the regions where 
the semiclassical approximation is not valid. 

A classical resonance in a perturbed integrable system manifests itself in the quantum 
picture if the size of the resonant island is greater than the Planck constant. In a billiard 
with a perturbed boundary this translates to the condition v LSL > A where L is the size of 
the billiard, SL ~ eL is the magnitude of the perturbation, and A = 1/fc is the wavelength. 
When this is the case, the resonant perturbation theory is required. The theory can be 
constructed by expanding the phase of the PSS wavefunction \p in ke M , starting with M = |. 
In the leading order the perturbed part of the wavefunction satisfies the one-dimensional 
Schrodingcr equation with the effective potential that scales as k 2 e. Because of the above 
condition, the equation has bound state solutions that indicate the localized two-dimensional 
wavefunction. Usually the localization is in the neighborhood of a stable periodic orbit. The 
two-dimensional wavefunction can be found by propagating tjj from the PSS into the bulk of the 
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system, basically along the classical paths. The spread of the wavefunction in the momentum 
space scales as k^/e. The non-resonant states are delocalized in the coordinate space, but are 
better localized (as fee) in the momentum space. Our theory automatically takes into account 
the whole perturbed family of periodic orbits, not just the neighborhood of a stable orbit. The 
solutions of the one-dimensional Schrodinger equation include both bound and unbound and 
thus cover the transition from the resonant to non-resonant states, as long as the influence 
of other resonances can be neglected. In this work we considered two examples of billiards 
with a perturbed boundary: the circle, in particular, the short stadium, and the rectangle, in 
particular, the tilted square. 

In some instances the perturbation theory can be applied to the non-perturbative systems. 
This can be done if the system is close to some integrable system in the regions where the 
wavefunction is large. We mentioned the ice-cream cone billiard that has the localized states in 
the circular part, the bouncing ball states in the stadium localized outside of the semicircles, the 
bouncing ball states near a period-2 orbit, and the whispering gallery modes. The system may 
have an additional perturbation. For example, the bouncing ball modes in the slanted stadium 
are shifted towards the wider end. The whispering gallery mode in a convex billiard follows 
from the standard EBK quantization procedure of the respective classical motion. However, 
when the boundary has a point or a region of zero curvature, there is no classical adiabatic 
invariant and the classical particle is not localized near the boundary. Still, the localized 
quantum state may exist. Outside of the region of zero curvature it is a superposition of 
the standard whispering gallery modes that are well-defined. The singular region acts like 
a scatterer that mixes the standard modes. The scattering matrix can be approximately 
expressed in terms of the T-operator. If the cigcnstatcs of the scattering propagator are 
composed of the small number of the standard modes, they are localized near the boundary. 

The magnetic field is a small perturbation for charged particle in a square billiard if its 
momentum p ^> eA/c, where A is the vector potential. Then, in the uniform field case, 
the cyclotron radius is greater than the size of the square. Thus, in the leading order, the 
shape of the orbits is unchanged. The resonant perturbation theory should be applied when 
kL^l = ^/{e/cfijkBL 3 > 1 where B is the magnetic field. In this case there exist states that 
are localized near the stable periodic orbits. The unstable periodic orbits are the time-reversals 
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of stable orbits. The states associated with the (1,1) resonance dominate the susceptibility. 
Since the time-reversal symmetry is broken, the states carry the persistent currents. The low 
transverse modes are paramagnetic and the higher modes are diamagnetic. In the square with 
a flux line the localization also takes place but the diffraction effects are strong. The diffraction 
is smaller for a finite size tube. 

When one of the system's coordinates changes much faster than the other, one can solve 
the two-dimensional Schrodingcr equation in the Born-Oppenheimer approximation (BOA). 
The solution usually agrees with the perturbation theory in the leading order. The standard 
example are the bouncing ball states where there is a natural separation of the fast and slow 
motion. The rectangular geometry is easy to deal with as well, if one employs the method of 
images. In this case the BOA is similar to the channeling approximation. On the other hand, 
in the circular geometry only the low angular momentum states are suitable for the BOA. 
Even then it fails in the polar coordinates near the center. Thus the BOA solution cannot be 
completed in this case since it lacks the boundary condition at the center. At the same time, 
it agrees with the perturbation theory in the asymptotic region if the missing information is 
provided. In the whispering gallery mode the slow variable is the curvature of the boundary. 
In this instance the BOA solution is valid also near the caustic and in the classically forbidden 
region where the perturbation theory is not applicable. 

The T-operator contains all the information about the energy spectrum in the semiclassical 
approximation. In particular, the oscillatory part of the density of states can be expressed in 
terms of the traces of the powers of T-operator. Depending on the dynamics of the system 
the traces can be evaluated to yield the Gutzwiller trace formula in the hard chaotic case 
or the Berry- Tabor formula in the integrable case. In the former the periodic orbits are 
isolated and in the latter they form families. The perturbed integrable system connects the 
opposite situations and provides an interpolation formula. The interpolation formula can be 
parametrized. The parametrization depends on four parameters: the perturbed part of the 
action and the monodromy determinant for the stable and unstable periodic orbits. When 
the action difference between the stable and unstable orbits is greater than h, these orbits are 
well isolated, when the actions are equal, the system is integrable. In the case of the coupled 
quartic oscillators the parametrization works very well. 



155 



The results of this work may be hard to verify experimentally, since it is the energy spec- 
trum, not the probability density, that is usually measured. Nevertheless, there are experi- 
mental techniques available. The simplest would be to measure the electric field distribution 
in a microwave cavity. The experiments have already been done for a chaotic shape cavity, 
so there is no principal difficulty to make it almost integrable. The electric field in a cavity 
and the electron wavefunction in a billiard satisfy the same differential equation. The mag- 
netic billiard can be modeled by adding a ferrite strip on the walls of the cavity. In this case 
the analogy holds only in the leading order in vector potential. It should also be possible to 
measure the electron wavefunction in quantum corrals, but the effect may be diminished due 
to a substantial leakage through the walls. Another class of experiments involves the surface 
waves in water, although here the dissipation may be a big problem. 
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